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ABSTRACT 


Procedures  are  outlined  for  solution  of  the  slowing-down  equation  for 
fast  electrons- by  direct  numerical  integration.  Spectra  for  H^  are  computed 
by  use  of  Gerhart's  cross  section  set,  suitably  modified,  and  with  inclusion 
of  vibrational  excitations.  The  flexibility  and  stability  of  the  approach  is 
demonstrated  by  solutions  for  source  energies  ranging  from  1 MeV  down  to  24  eV 
with  the  spectrum  extended  below  2 eV.  Energy  and  number  conservation 
are  checked  in  detail,  and  the  ionization  yield  is  computed  as  a function  of 
energy.  Spectrum  and  ionization  yield  so  obtained  are  compared  with  previous 
calculations;  general  agreement  is  found  with  some  significant  differences. 
Agreement  with  Combecher's  measured  ionization  yields  at  low  energies  is 
rather  satisfactory  overall,  but  with  significant  deviations  of  trend.  Some 
exploration  of  secondary  and  tertiary  electron  contributions  is  included. 
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NUMERICAL  SOLUTION  OF  THE  ELECTRON  SLOWING-DOWN  EQUATION  FOR  MOLECULAR  HYDROGEN* 

L.  V.  Spencer  and  R.  Pal 

I . Introduction 

Electron  energy  degradation  or  "si owl ng-down"  is  due  to  inelastic  inter- 
actions with  the  surrounding  medium.  The  process  generates  a spectrum  y 
which  expresses  the  relative  likelihood  that  a fast  electron  is  detected  in 
different  energy  intervals  dT.  More  specifically,  ydT  can  be  considered  as 
the  expected  value  of  the  pathlength  traveled  by  a source  electron,  together 
with  that  by  secondary  electrons  which  it  generates  during  the  slowing  down 
process,  while  the  one  or  any  of  the  others  have  kinetic  energy  within  dT. 

The  importance  of  y derives  largely  from  the  fact  that  one  need  only 
multiply  y_  by  the  (suitably  normalized)  cross  section  for  some  process,  and 
integrate  over  all  electron  kinetic  energies,  to  obtain  the  yield  for  that 
process,  whether  excitation  or  ionization.  But  given  a correct  spectrum  y, 
the  resulting  yields  are  no  more  accurate  than  the  cross  sections  which  are 
used  to  evaluate  such  integrals.  Errors  in  evaluation  of  the  slowing-down  can 
accentuate  and  greatly  complicate  the  errors  and  uncertainty  of  any  yields  so 
calculated. 

There  are  two  difficulties  associated  with  determining  the  energy  degra- 
dation spectrum.  First,  the  inelastic  scattering  cross  sections  are  complex 
and  hard  to  obtain.  Secondly,  because  of  the  complexity  of  the  cross  sections 
the  numerical  solution  of  the  transport  equation  is  not  a trivial  task.  As  a 
result  there  have  been  nemerous  investigations  with  various  types  of  approxi- 
mations which  will  be  briefly  reviewed  below.  But  since  the  basic  integral 
equation  of  which  y is  the  solution  can  in  principle  be  solved  by  standard 
numerical  techniques,  we  feel  that  it  is  appropriate  to  make  an 
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attempt  to  do  this.  This  paper  describes  a successful  calculation  of  y for 
molecular  hydrogen  by  application  of  well-known  numerical  integration  techniques 
to  the  standard  form  of  the  slowing-down  equation.  Because  the  cross  sections 
utilized  are  current  the  resulting  data  is  of  general  interest. 

Before  proceeding  to  describe  this  effort,  we  identify  some  of  the  other 
methods  which  have  been  applied  in  the  following  section,  two  of  the  most 
recent  (for  H2)  being  due  to  Douthat  [1]  and  to  Garvey,  Porter,  and  Green  [2]. 
Our  comparisons  with  results  from  these  and  other  calculations  are  of  interest 
perhaps  more  because  different  cross  section  sets  have  been  applied  than 
because  of  contrasting  computational  methods.  We  use  the  cross  section 
compilation  of  D.  E.  Gerhart  [3];  but  we  must  modify  Gerhart's  forms  to  render 
them  more  readily  applicable  to  our  procedures,  as  well  as  to  make  extensions 
to  relativistic  energies.  In  addition,  we  take  account  of  vibrational  modes 
in  order  to  extend  our  results  well  below  the  threshold  for  electronic 
excitations,  i.e.  down  to  about  1.5  eV. 

We  believe  that  the  calculations  described  here  are  significantly  more 
accurate  than  previous  evaluations  of  y for  molecular  hydrogen.  Comparison 
with  recent  measurements  of  the  ionization  yield  as  a function  of  electron 
source  energy  can  be  construed  as  pointing  out  the  importance  of  more  accurate 
ionization  and  excitation  cross  sections  in  the  range  from  15  eV  to  perhaps 
30  eV  for  future  improvement  of  both  slowing-down  spectrum  and  yield  results. 

New  aspects  explored  here  include  careful  and  complete  checks  of  energy 
conservation  and  conservation  of  the  number  of  primary  particles.  In  addition, 
the  spectral  structure  near  the  source  energy  is  resolved  in  its  most 
important  features;  and  there  is  some  discussion  of  general  character! sties  of 
spectral  structure,  particularly  as  it  results  from  those  generations  of 
secondary  electrons  which  are  its  most  important  components. 
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1 1 . Comments  on  Computational  Methods 


If  electrons  are  generated  with  initial  kinetic  energy  T^,  the  slowing- 
down  equation  of  which  y is  the  solution  is 

K(T)  y(T,Tj  = ^ dr  k(T',T)  y(TM^)  + 6(T-TJ  , (1) 

where  K(T)  is  the  total  inelastic  interaction  probability  per  unit  pathlength 
at  kinetic  energy  T;  and  k(T',T)dT  is  the  differential  probability,  per 
unit  path  length,  for  an  electron  of  energy  T'  to  undergo  an  interaction 
resulting  in  an  electron  having  kinetic  energy  in  the  interval  dT  containing 
T.  This  electron  can  be  a new  secondary,  or  the  result  of  an  energy  loss 
(T'-T).  The  Dirac  delta  function  in  eq  (1)  describes  a monoenergetic  source 
of  unit  strength  with  kinetic  energy  T^;  and  the  spectral  function  y(T,T^) 
has  dimensions  path  length  per  unit  energy. 

Direct  numerical  integration  of  eq  (1)  has  not  been  thought  very  feasible 
because  of  the  general  structure  of  k(T',T)  for  T'  near  T.  A large  and 
complex  concentration  of  interaction  probabilities  occupies  a narrow  region 
11  eV  < T'-T  < 18  eV,  as  seen  in  Fig.  1.  Further,  for  whatever  value  of  i 
may  be  selected,  energies  within  a few  tens  of  eV  of  T^  can  only  be  reached 
if  no  more  than  one,  two,  or  at  most  a few  interactions  have  occurred.  This 
means  that  the  structure  for  T'  near  T in  the  cross  sections  also  appears 
for  T near  T^  in  y(T,T^).  Thus  the  integrand  of  eq  (1)  has  complicated 
features  for  T'  near  T,  and  for  T'  near  T^;  and  these  two  regions  of 
structure  both  make  significant  contributions  to  the  integral.  Direct 
numerical  integration  must  therefore  deal  with  the  necessity  to  use  integration 
steps  perhaps  less  than  1 eV  wide  to  resolve  structure,  while  at  the  same 
time  covering  an  energy  range  below  T^  that  can  be  of  the  order  of  10^  eV 
or  greater. 
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Most  calculations  with  above  a few  hundred  eV  have  attempted  to 
bypass  the  complexities  of  y(T5T^)  for  T near  T^.  Incidentally,  this 
structure  was  first  described  by  H.  W.  Lewis,  and  therefore  has  been  called 
the  "Lewis  effect"  [4].  The  various  computational  methods  fall  into  distinct 
classes  as  follows; 

(1)  Numerical  evaluation  of  an  integral  form  of  eg  (1).  If  one 

integrates  all  terms  of  eq  (1)  over  the  variable  T,  from  a value  (say)  , 
to  the  source  energy  ‘’'q,  one  obtains  an  equation  of  balance  for  the  number 
0-^  electrons  with  kinetic  energy  above  T.  ^ A similar  integral  over  al'' 
terms.,  in  which  an  additional  weighting  factor  T is  used  in  the  resulting 
integrands,  gives  a corresponding  equation  of  balance  or  "conservation"  for 
energy  flowing  into,  and  out  of,  the  reservoir  of  electrons  with  kinetic 
energies  between  T and  These  conservation  equations  exhibit  the 

significant  advantage  of  having  kernel  functions  with  less  irregular  and 
singular  trends  than  eq  (1).  Further,  conservation  of  number  or  energy  is 
usually  more  important  than  a detailed  description  of  cross  section  structure, 
when  one  attempts  to  solve  for  y at  most  spectral  energies.  A number  of 
authors  have  successfully  used  the  "number  balance"  equation  (eq  (7))  to 
obtain  solutions  for  y [1,5-10]. 

(2)  A more  elementary  and  older  method,  which  is  sometimes  used  to 
obtain  approximate  solutions  to  the  slowing  down  equation  employs  the  con- 
t j nuous  slowing  down  approximation  (CSDA)  [11 -1 5].  While  the  CSDA  is  more 
suited  to  the  slowing  down  of  heavy  charged  particles,  it  also  describes  in  a 
rough  way  the  consequence  of  the  frequent,  small-loss  interactions  of  electrons 
with  the  medium  through  which  the  electrons  move.  One  can  readily  derive  this 

’■See  eqs  (7)  and  (8). 
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approximation  from  eq  (1).  Any  electron,  once  generated,  contributes  a term 
1/S(T)  to  the  slowing-down  spectrum  below  the  energy  at  which  generated, 
where  S(T)  is  the  stopping  power;  and  the  number  of  secondary  electrons 
generated  can  be  calculated  from  a much  simplified  integral  equation.  Calcu- 
lations of  this  type  give  a good  general  description  of  y(T,T^)  which  can  be 
sufficiently  accurate  for  some  purposes.  They  are  easily  developed  and 
inexpensive,  but  they  are  not  highly  accurate. 

(3)  Monte  Carlo  methods  can  be  used  to  evaluate  the  slowing-down  spectrum 
as  an  average  over  a few,  or  perhaps  a few  hundred  or  thousand  case  histories 
[16-19].  All  that  is  required  for  good  statistics  is  many  collisions  in  each 
energy  interval  specified.  Since  each  source  electron  not  only  remains  a part 
of  the  calculation  to  very  low  energies,  but  also  gives  rise  to  many  sec- 
ondaries, the  statistics  actually  improve  as  the  calculation  progresses  to  low 
energies.  Energy  conservation  can  be  readily  imposed  by  insuring  that  all 
energy  resulting  from  each  interaction  is  properly  accounted  for.  The  Lewis 
effect  region  can  be  described  if  and  as  desired.  Hence  results  obtained  by 
use  of  this  method  can  be  impressive.  Perhaps  its  main  drawback  is  the  time 
required  to  calculate  each  case  history,  which  makes  computations  of  many  case 
histories  with  high  source  energies  comparatively  expensive. 

(4)  Use  of  a discrete  analog,  in  which  electrons  progress  from  higher  to 
lower  energy  "bins"  according  to  computed  transition  probabilities  between 
pairs  of  bins.  Mostly,  such  calculations  have  been  used  with  rather  low 
source  energies  (e.g.,  < 1 keV)  [20-24].  A recent  publication  gives  results 
of  systematic  calculations  of  this  type  to  much  higher  energies  [2].  The 
authors  of  this  last  paper  chose  an  interval  near  the  source  energy  suffi- 
ciently wide  that  the  Lewis  effect  was  wiped  out;  but  they  indicate  that  the 
method  can  be  extended  to  give  a more  detailed  description  if  desired. 
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(5)  Use  of  methods  which  solve  an  equation  adjoint  to  eq  (1)  [3,24-28]. 
Such  equations  have  a form  similar  to  eq  (1);  but  the  solution  progresses 
upward  in  energy  rather  than  downward.  Difficulties  with  the  integration 
process  are  not  removed,  because  the  complexities  of  the  function  k(T',T) 
remain.  These  calculations  usually  evaluate  yields  directly,  rather  than 
y(T,T^).  Most  such  calculations  have  been  performed  at  low  energies  (up  to 

1 keV);  and  the  development  of  a wide  variety  of  methods  to  solve  these 
equations  is  still  in  process. 

(6)  The  application  of  numerical  integration  to  eq  (1),  which  is  the 
subject  of  this  paper,  and  which  has  greatest  resemblance  to  method  (4),  that 
is,  use  of  a discrete  analog.  Here  one  attempts  to  apply  standard  techniques 
of  numerical  integration  to  the  evaluation  of  the  scattering  integral,  for  T 
values  which  progress  stepwise  to  lower  energies,  some  use  of  this  method  at 
low  energies  is  reported  by  Douthat  [1]. 

We  know  of  two  attempts  to  make  comparisons  of  the  different  methods  when 
the  same,  or  closely  similar,  cross  sections  are  being  used  [23,24].  The  more 
recent  study  features  W,  the  reciprocal  of  the  ionization  yield,  as  an 
integral  comparison  parameter  [24].  It  gives  little  information  on  relative 
accuracies  of  different  types  of  computations  of  y(T).  The  earlier  study 
[23]  compared  CSDA  results  for  He  and  N2  with  results  obtained  by  the 
discrete  analog  method.  Only  electron  sources  of  500  eV  and  below  were  used. 
Differences  are  significant;  but  qualitative  and  rough  quantitative  agreement 
are  found.  One  reason  for  attempting  to  apply  direct  numerical  integration 
systematically  is  the  possibility  that  accuracy  for  a given  cross  section  set 
can  thereby  approach  ~ 1%.  This  would  effectively  eliminate  computational 
methods  as  a source  of  significant  error  in  evaluation  of  si owi ng-de-wn  spectra 
or  yields.  Such  accuracies  probably  require  a more  detailed  evaluation  of  the 
Lewis  effect  region;  but  this  can  be  desirable  on  its  own  merits. 
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A certain  amount  of  structure  occurs  at  low  energies  in  the  si owi ng-down 
spectrum  also,  because  the  thresholds  for  different  types  of  transition  can 
give  rise  to  irregularities.  We  think  that  this  structure  can  be  more  easily 
evaluated  by  direct  integration  of  eq  (1)  than  by  use,  e.g.,  of  integral  forms 
of  the  transport  equation  (method  (1)  preceding). 

Successful  methods  for  direct  numerical  integration  also  can  in  principle 
be  used  to  evaluate  space  and  directional  aspects  of  electron  transport  by 
methods  not  heretofore  applied  with  realistic  cross  section  data  to  this  type 
of  particle  transport. 

Because  the  systematics  of  slowing-down  spectra  are  important  for  obtaining 
similar  spectral  information  about  materials  and  source  energies  not  otherwise 
studied  in  the  same  detail,  we  have  organized  the  calculations  to  permit 
identification  of  the  primary  electrons  as  well  as  the  1st,  2nd,  and  higher 
generations  of  secondary  electrons.  The  different  generations  play  rather 
different  roles;  and  their  spectral  components  require  a somewhat  different 
description  for  each  generation. 

One  problem  avoided  in  our  calculations  for  ^2  is  the  problem  of  inner 
shells.  While  hydrogen  has  only  one  shell,  most  other  materials  of  interest 
have  inner  shells.  In  the  kernel  function  k(T',T)  different  atomic  or 
molecular  shells  can  contribute  terms  describing  electrons  with  the  same 
kinetic  energy  T;  so  one  must  in  principle  include  calculations  for  all 
contributing  shells  whenever  k is  required.  These  calculations  must  be 
performed  many  times,  because  k appears  as  a factor  at  each  step  of  each 
integration.  We  do  not  see  any  fundamental  difficulty  in  this  requirement, 
although  it  surely  increases  the  computational  time  of  calculations  by  any 
method. 
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A related  problem  which  likewise  appears  only  in  materials  with  Z > 2 
is  that  of  including  Auger  electron  components.  These  have  energies  character- 
istic of  transitions  between  different  shells.  They  must  be  treated  as  an 
additional  source  whose  inclusion  is  necessary  for  conservation  of  energy. 
Again,  there  appears  to  be  no  fundamental  difficulty  in  accounting  for  this 
component;  but  our  current  task  is  simplified  by  their  absence  from 
spectra. 


III.  Cross  Sections  and  Stopping  Powers 

The  shape  of  the  differential  cross  sections  has  a dominant  influence  on 
integration  procedures  appropriate  for  this  and  other  similar  problems. 

Section  IV,  which  follows,  illustrates  this  in  detail.  In  this  section  we 

give  general  information  about  the  set  of  differential  cross  sections  with 
which  we  work,  and  the  total  cross  section  and  stopping  power  data  derived 
from  it.  More  details  about  these  cross  sections  are  given  in  Appendix  I. 

We  have  adopted  a well-known  set  of  cross  sections  for  our  calculations, 
namely  the  cross  section  set  developed  by  Donald  E.  Gerhart  in  his  PhD 
thesis  [3].  Gerhart  followed  Platzman's  use  of  the  Weiszacker-Wil 1 iams 
approach  in  which  collisions  are  classified  as  either  hard  or  soft,  but  not 

both.  This  involves  specifying  cut-offs  to  both  soft  ("glancing")  and  hard 

("knock-on")  collisions. 

In  the  process  of  adapting  these  cross  sections  to  our  use,  we  have  made 
a number  of  modifications.  A brief  summary  of  these  changes  is  as  follows: 

1)  We  have  extended  Gerhart's  cross  sections  into  the  relativistic  range  in  a 
manner  which  gives  stopping  powers  of  the  Bethe  theory  type,  with  reasonable 


values  for  the  "mean  excitation  potential",  here  designated  (I).  2)  In 

addition,  we  have  modified  the  form  of  Gerhart's  cut-off  function  for  small 
energy  losses  due  to  hard  collisions;  for  our  purpose  it  is  more  advantageous 
to  use  an  abrupt,  rather  than  gradual  cut-off  of  this  kind.  3)  At  low  electron 
kinetic  energies,  where  the  differential  spectrum  of  energy  losses  is  poorly 
known, extend  the  use  high  energy  cross  section  forms  by  ad  hoc  adjustments 
of  the  cut-off  for  momentum  transfers,  which  is  a constant  at  high  energies 
in  cross  section  formulae  for  glancing  collisions.  This  replaces  an  abrupt 
transition  at  350  eV  to  low  energy  forms  which  Gerhart  used,  while  preserving 
good  agreement  with  total  ionization  cross  sections.  4)  We  have  added  vibra- 
tional excitation  cross  sections  which  were  omitted  by  Gerhart,  and  which  play 
a dominant  role  for  kinetic  energies  below  about  12  eV.  5)  Finally,  we 
symmetrized  the  glancing  collision  cross  sections.  These  matters  are  given 
in  greater  detail  in  Appendix  I,  which  also  contains  additional  numerical 
information  omitted  from  Gerhart's  publication,  such  as  a set  of  transition 
energies  for  the  discrete  part  of  the  spectrum  as  given  by  Namioka  [29]. 

In  the  process  of  satisfying  sum  rules,  Gerhart  adopted  an  oscillator 
strength  distribution  (OSD)  which  gives  a value  for  the  ionization  dipole 
strength  ~ 1 OX  lower  than  data  from  other  sources,  with  OSD  values  for  the 
excitation  region  correspondingly  a little  higher  than  other  data.  We  doubt 
that  this  is  desirable,  and  one  of  us  (RP)  has  developed  an  OSD  safisfying 
the  sum  rules  without  such  adjustments.  But  the  importance  of  Gerhart's  cross 
sections  is  such  that  we  have  used  them  --  with  modifications  as  outlined 
above  — in  the  extensive  numerical  explorations  reported  here. 
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Figures  1 and  2 present  differential  cross  section  data  as  calculated  by 
our  version  of  Gerhart's  cross  sections.  Figure  1 displays  the  variation  with 
energy  loss.  Figure  2 shows  that  our  version  of  Gerhart's  cross  sections 
agrees  rather  well  in  shape  with  the  Opal,  Beaty,  and  Peterson  [30]  experi- 
mental data  for  the  spectrum  of  secondary  electrons  at  T = 500  eV.  One 
should  note,  however,  that  the  experimental  data  have  been  divided  by  a scale 
factor  of  1.65.^  A similar  comparison  at  T = 100  eV  with  recent  experimental 
data  due  to  Rudd  and  Dubois  [31]  shows  somewhat  poorer  agreement  as  to  shape, 
but  better  agreement  in  magnitude.  In  both  cases,  the  experimental  data 
appear  to  drop  below  the  theoretical  data  at  very  low  secondary  electron 
energies. 

Figures  3 and  4 present  stopping  power  and  total  cross  section  results. 
Values  for  the  stopping  power  obtained  using  Bethe's  expression  are  also  given 
in  fig.  3,  as  are  some  recent  unpublished  data  due  to  Green  [32].  Gerhart's 
OSD  leads  to  a computed  value  for  <I)  of  about  19.2  eV  as  stated  in  his 
paper  [3];  and  this  value  was  used  in  calculating  the  curve  labeled  "Bethe". 

On  the  other  hand,  our  set  of  cross  sections  leads  to  an  equivalent  constant 
with  the  value  <I)p  = 18.6  eV  even  though  utilizing  Gerhart's  OSD.  While  it 
is  possible  to  make  adjustments  which  would  increase  the  latter  value,  we  have 
not  attempted  to  do  so.  See  Table  1.6  for  other  data  on  <I>,  both 
experimental  and  theoretical. 


^This  means  that  the  Opal,  Beaty,  Peterson  data  extrapolated  in  a reasonable 
way  give  an  integrated  value  significantly  larger  than  data  due  to  Rapp  and 
Englander-Golden  [33],  to  which  Gerhart's  cross  sections  are  normalized. 
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IV.  General  Method  of  Solution 


A.  Comments  on  the  integration 

In  solving  integral  equations  such  as  eq  (1),  the  usual  approach  involves 
replacement  of  the  integral  by  a weighted  sum  as  follows: 

dT  k(T',TJ  y(T',Tjj)  - ^ w-  k(T.,T^)  y(T.,T^)  . (2) 

If  all  values  y(T.,T^)  are  known  except  for  y(T  ,T  ),  and  if  all  values 
for  w^.  have  been  specified,  eq  (1)  gives  a value  for 

When  both  factors  of  the  integrand  vary  smoothly  and  slowly  from  one 
value  of  T^.  to  its  neighboring  values,  weights  w.  can  be  adopted  from 
one  of  the  standard  formulae,  such  as  Simpson's  rule  or  even  the  trapezoidal 
rule.  This  method  of  calculation  is  suitable  for  regions  in  which  k(T',T) 
and  y(T,T^)  vary  with  similar  rapidity  and  irregulari ty , because  the  many 
values  for  X i^equired  to  describe  such  functions  correspond  rather  well  in 
spacing  and  number  to  the  mesh  of  values  required  to  evaluate  the  scattering 
integral  accurately. 

An  altogether  different  circumstance  can  occur,  however.  One  of  the 
factors  in  the  integrand  can  be  rapidly  varying  while  the  other  is  smooth, 
slowly  varying,  and  predictable.  Thus  for  large  T^,  at  most  kinetic  energies 
T',  T'  < T^,  one  has  a slow  and  predictable  function  y(T',T^)  combined 
with  very  rapid  and  irregular  variations  of  k(T',T^)  whenever  T'  is  within 
perhaps  30  eV  of  T^.  This  is  wel 1 -i 1 1 ustrated  in  fig.  1,  where  the  cross 
section,  which  is  proportional  to  k,  is  given  as  a function  of  the  difference 
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between  the  variables.  The  irregular  behavior  for  (T^  - T)  or  (T'  - T^) 
between  11  eV  and  25  eV  occurs  for  virtually  all  T^,  regardless  of 
trends  of  y. 

Under  these  circumstances  eg  (2)  is  still  appropriate;  but  it  is  no 
longer  desirable  to  use  elementary  numerical  integration  weights.  Instead, 
one  can  use  integrals  over  k,  usually  moments,  to  construct  an  approximation 
to  k as  a sum  of  Dirac  delta  functions.  The  integral  is  then  readily 
performed  and  gives  an  expression  of  the  type  shown  by  eq  (2).  We  refer  here 
to  such  "6  - function  approximations"  by  the  general  designation  "Gauss 
quadrature"  (GQ) . ^ 

The  €imple  numerical  integration  algorithms  for  the  solution  of  eq  (1)  on 
the  basis  of  expressions  of  the  type  illustrated  by  eq  (2)  require  the  T^  to 
have  the  same  values  for  different  T^.  Other  algorithms  are  possible;  but 
this  is  preferable  on  the  basis  of  simplicity.  We  thus  choose  T^.  values  for 
the  complete  calculation  at  the  start,  with  T^ ' s which  march  down  this 
preselected  list.  To  achieve  accuracy  comparable  with  Simpson's  rule  in 
eq  (2),  we  evaluate  three  energy  moments  of  k for  the  combined  pair  of 
intervals  T^_^  < T'  < T^.  The  resulting  three  conditions^  fix  6-function 
weights  at  T^,  , T^  These  weights  can  be  cumulated  with  other 

weights,  since  all  play  the  same  role  in  the  numerical  integration  at  various 
T^.  values  of  the  same  set. 


^We  here  include  the  option  of  pre-specifying  some  or  all  x.  values.  This 
results  in  a unified  view  of  numerical  integration  procedures  of  many  types. 
The  usual  Gauss  quadrature  can  be  considered  as  an  approximation  of  the  type 

N 

w(x)  ==  a^6(x-x^.  )>  where  the  a^.  and  x^  are  fixed  by  2N  jrioment 

integrals.  In  the  classic  case,  the  weight  function  w(x)  is  unity  for 
-1  < X < 1.  For  another  well-known  case  w(x)  = e~^,  with  0 < x < 
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By  restricting  our  approach  to  integration  formulae  which  yield  sums  of 
the  type  shown  in  eq  (2),  we  can  determine  w_j  values  by  two  or  more  different 
procedures  of  this  type,  each  appropriate  to  a sub-region  of  the  domain  of 
i ntegration. 

B.  The  integration  mesh 

Description  of  the  Lewis  effect,  i.e.  resolution  of  complex  spectral 
details  near  the  source  energy,  requires  that  intervals  for  numerical 
integration  be  ^ 1 eV,  as  already  noted  (see  fig.  1).  Similarly,  at  low 
energies,  near  and  below  the  ionization  potential,  trends  of  y(T,T^)  with  T 
are  again  steep,  with  possible  structure.  Small  intervals  are  therefore  also 
appropriate  for  solution  of  the  integral  equation  at  low  energies.  But  it  is 
quite  clear  that  most  of  the  energy  range,  that  between,  say,  100  eV  and 
(T^  - 100)  eV,  does  not  require  such  a detailed  description  because  y(T,T^) 
has  trends  which  are  smooth  and  not  particularly  rapidly  varying. 

Our  choice  of  the  energies  constituting  the  solution  mesh  was 

therefore  as  follows:  We  used  the  tabulation  intervals  established  for  the 

OSD  to  develop  a set  of  values  in  the  energy  range  from  down  a total 

of  40  to  80  eV  lower  in  energy.  The  interval  size  was  then  increased  abruptly 
from  AT  ~ 2 eV  to  AT  = 20  eV.  Following  this,  the  intervals  were  increased 
in  geometric  progression  until  the  energy  loss  amounted  to  T^/2.  Typically, 
then,  the  energies  T^,  such  that  T^  < T^/2,  were  equated  to  the  energy 
1 osses  (T^  - T^).  Thus  the  solution  energies  were  approximately  twice  in 
number  those  required  to  carry  the  solution  from  T^  to  T^/2. 

The  reason  for  an  abrupt  change  in  interval  size  as  referred  to  above  is 
that  intervals  smaller  than  2 eV  can  resolve  at  least  moderately  well  the 
main  features  of  the  glancing  collision  region,  even  between  13  and  17  eV. 

On  the  other  hand,  intervals  greater  than  20  eV  are  large  enough  so  that  the 
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whole  band  of  irregularities  in  k(T',T)  is  contained  within  a single 
interval;  and  fine  resolution  is  then  rendered  unnecessary  except  in  regards 
to  the  calculation  of  GQ  weights,  which  is  accomplished  by  interpolation  on 
basic  tables  of  cumulative  integrals. 

We  did  find  it  convenient  to  buffer  the  transition  to  large  intervals 

with  two  intervals  of  size  10  eV,  at  the  same  time  ensuring  that  the  GQ 

region  extended  beyond  one  such  interval.  Table  1 gives  an  illustrative  set  of 

T values,  for  a 1 keV  source, 
n 

From  the  attention  given  to  such  details,  it  should  be  evident  that 
setting  up  a suitable  integration  mesh  constitutes  a significant  part  of  our 
calculation.  Typically,  the  fine  spacing  for  high  and  low  energy  regions 
requires  about  40-45  points  each,  with  up  to  90  additional  points  assigned 
to  the  remaining  (usually  large)  energy  region.  We  used  up  to  ~ 200 
integration  points,  so  that  instruction  and  data  memory  components  together 
occupied  just  under  64  K cells  of  the  memory. 

C.  The  integration  procedure 

Both  glancing  and  knock-on  cross  section  components  have  a rather  simple 
analytic  structure;  thus  k(T',T),  with  T'  = T + E,  has  the  form 

k(T+E,T)  = E fyi+E)  gyi.E)  . (3) 

i 

The  functions  f^.  (T  + E)  are  slowly  varying,  and  the  functions  g^-(T,E)  can 
either  be  integrated  analytically,  or  depend  only  on  E,  or  both.  We  have 
already  mentioned  tabulations  of  cumulative  integrals  of  the  (rapidly  varying) 
g^.(E)  for  glancing  collisions.  Corresponding  functions  g^-(T,E)  for  knock- 
on  collisions  can  be  integrated  analytically. 
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Because  the  functions  E)  are  slowly  varying,  if  combined  with 

the  function  y(T  + E,T^)  the  product  tends  to  vary  slowly  or  rapidly  according 
to  the  characteri sties  of  y.  This  leads  us  to  pretabulate  the  f^.  at  the 
solution  energies  T^,  with  integration  based  on  variables  G^.  defined  by 

dG 

3r=9i(T>E)  • (9) 


Then 


k(T',T)  y(T',iy  dT' 


E JdG,  f,(T,G, 


1 1 


y(T,G.,T^) 


(4') 


This  approach  should  be  contrasted  with  the  more  direct  procedure  of 
simply  evaluating  k(T',T)  as  necessary  for  all  combinations  (Tj,T^) 
which  appear  in  the  set  of  integrations.  We  have  ordinarily  used  5 or  6 
terms,  i.e.  values  of  i,  so  that  the  right  side  of  eq  (4')  has  ~ 6 terms 
in  the  sum.  On  the  other  hand,  pre-tabulation  of  the  f^.  is  quick  and  cheap, 
because  it  need  be  done  only  once  to  obtain  a list  of  values  for  all  solution 
energies  T..  By  contrast,  direct  evaluation  of  k(T',T)  requires  calcula- 
tions  of  all  component  factors  of  k for  a quadratic,  rather  than  linear 
array,  more  or  less  cancelling  out  any  extra  computer  time  required  by  the 
multiplicity  of  i values.  In  addition,  by  differencing  G^.  we  readily 
obtain  the  lowest  moment  required  for  the  GQ;  repetitious  tabulations  to 
evaluate  GQ  moments  are  avoided  by  use  of  several  such  cumulative  functions. 
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We  believe  that  it  would  be  about  equally  feasible  to  use  tabulations  of 
cumulative  integrals  for  the  Gauss  quadrature,  with  direct  evaluation  of 

k(T',T)  reserved  for  evaluation  of  the  remainder  of  the  scattering  integral; 
but  we  have  thus  far  chosen  not  to  proceed  in  this  way. 

We  always  use  GQ  for  the  final  interval,  T^_.|  > T'  > T^,  and  we  use  a 
version  of  Simpson's  rule  for  the  remainder  of  the  scattering  integral.  Our 
GQ  contribution  to  eq  (2)  is  as  follows:‘" 

T. 

dT'  k(TM„)  y(TVy  =»  W,y(Tj^,,Ty  + -^0) ' ( = ) 

n 

where  j_  is  commonly,  but  not  always,  equal  to  n-1 . Since  k(T',T)  vanishes 
for  (f  - T)  < 11.2  eV,  we  find  it  advantageous  to  extend  the  GQ  interval 
to  fully  cover  this  gap,  because  other  integration  procedures  are  rendered  • 
inaccurate  when  extended  into  a region  in  which  the  integrand  is  zero.  The 
rule  adopted  for  choice  of  ^ is  therefore  that  it  equal  the  largest  m 
value  such  that  (T^  - T^)  > 11.2  eV. 

Six  cumulative  integrals  over  the  OSD  giving  data  for  three  energy 
moments  (see  Appendix  I),  together  with  corresponding  analytic  integrals  over 
components  of  the  knock-on  cross  section,  give  moments  of  the  different 
g^.  (T,E)  from  which  integration  weights  can  be  systematically  developed  for 
all  j_  values  by  combining  contributions  from  adjacent  integration  intervals, 
each  covering  two  intervals  of  the  tabulation. 


“^Note  that  T.  need  not  be  included  between  T and  T-. 

J- 1 n j 


/ 
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We  should  comment  that  in  integrals  of  the  type  shown  on  the  right  of 
eq  (4‘),  for  all  ranges  except  those  where  6Q  is  applied  we  systematically 
use 


dG. 


dG. 



du 


du  , 


(6) 


where  for  each  three  consecutive  T.j  values  the  corresponding  u_  values  are 
-1,  0,  and  1,  thus  permitting  application  of  the  elementary  form  of 

Simpson's  rule.  We  determine  dG^./du  by  fitting  each  three  consecutive  G^. 
values  to  a quadratic  in  u.  Note  that  only  two  of  the  six  cumulative  integrals 
over  the  g^-  are  required  for  this  process;  we  use  the  other  information 
available  in  Table  1.2  only  for  the  GQ. 

D.  Energy  and  number  conservation 

In  order  to  monitor  the  calculations  very  closely  we  wrote  a special 
subroutine  (ECHEK)  to  perform  energy  conservation  checks  at  selected  energies 
(see  Appendix  II).  The  integral  forms  of  the  slowing-down  equation  which 
were  referred  to  in  the  introductory  discussion  are  (for  energy  conservation) 

Tq 

dT'  S(T',T)  y(T',Tj  = Tq  , (7) 

T 

where  S(T',T)  is  a limited  stopping  power  which  excludes  energy  transfers  to 
secondary  electrons  with  kinetic  energies  exceeding  T.  For  number  conserva- 
tion the  equation 


2T+I  T^ 

J dT'  Kp(T',T)  y(TM^)  = 1 + r dT'  K^T’ ,T)  y(T' ,T^)  (8) 

' 2T+I 
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holds,  where  Kp(T',T)  and  K^(T',T)  are,  respectively , the  probabilities 
per  unit  path  length  for  removal  and  addition  of  an  electron  to  the  kinetic 
energy  interval  between  T and  T'  [5].  Accurate  evaluation  of  the  integrals 
in  eqs  (7)  and  (8)  is  of  comparable  difficulty  to  the  main  problem  of  evalu- 
ating eq  (1),  although  one  can  use  devices  that  have  previously  proven 
satisfactory  [5-10].  The  accuracy  of  the  numerical  integration  remains  a 

problem  because  S(T',T)  and  K (T',T)  give  important  contributions  for 

P 

(T'  - T)  small,  while  retaining  some  irregularity. 

For  evaluation  of  eqs  (7)  and  (8)  we  wrote  a general  stopping  power 
subroutine  (STN)  and  an  integrated  cross  section  subroutine  (TCX)  according  to 
a pattern  which  permitted  rapid  and  simple  cumulative  integration  over 
transition  energies  up  to  an  arbitrary  limit  (see  Appendix  II).  These 
routines  were  essential  to  the  evaluation  of  the  kernel  functions  of  eqs  (7) 
and  (8).  In  preparing  the  stopping  power  subroutine,  all  contributions  due 
to  the  energy  of  ionization  must  be  sometimes,  though  not  always,  included. 

This  requires  a special  term  when  fast  (i.e.  above  cut-off)  secondary 
electrons  were  generated. 

In  the  use  of  both  stopping  power  and  integrated  cross  section  sub- 
routines, one  must  bear  in  mind  that  for  T'  < 2T^  + I both  primary  and 
secondary  electrons  from  an  interaction  can  lie  below  the  energy  cut-off  T^. 
For  such  interactions  the  total  energy  of  the  electron  being  scattered  is 
"dissipated. " 

The  exercise  of  taking  account  of  energy  conservation  has  proved^  necessary 
and  very  instructive.  By  this  means  we  found  that  in  the  transition  from 
finely  spaced  integration  intervals  to  coarse  integration  intervals,  some  of 
the  energy  (typically  ~ 4%)  can  be  lost.  This  may  be  partly  due  -to  the  fact 
that  the  changeover  necessarily  occurs  sufficiently  near  the  source  energy  to 
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be  affected  by  the  Lewis  effect  irregularities.  The  first  few  points  with 
coarse  spacing  are  uncertain  by  a few  percent  for  this  reason;  and  any  resulting 
error  is  then  "locked"  in,  because  the  major  contribution  to  the  scattering 
integral  is  always  due  to  a few  points  thus  affected. 

To  correct  such  errors,  we  have  calculated  energy  conservation  as  well  as 
a value  for  y^.  at  the  first  five  T^.  values  with  the  coarse  spacing.  We 
adjusted  these  y.  values  immediately  as  if  any  deficiency  of  energy  conser- 
vation were  solely  due  to  the  y^-  value  being  determined.  This  worked  well. 

In  general,  we  find  that  energy-conservation  results  give  a good  check  on 
the  accuracy  of  the  computations.  The  number  conservation  data  are  in 
numerical  agreement  with  energy  conservation  data  above  T^/2.  Below  T^/2 
the  number  conservation  data  simply  record  the  total  number  of  fast  electrons 
generated  by  all  processes.  But  the  number  conservation  check  can  be  extended 
to  very  low  energies  by  applying  it  to  the  primary  electron  component  only. 

E.  Special  features 

At  the  time  the  computer  program  was  being  written,  we  felt  that  it  might 
be  very  important  to  distinguish  between  primary  electrons,  the  first  generation 
of  secondary  electrons,  and  higher  generations  of  secondaries.  Hence  in 
setting  up  the  basic  table  of  integration  weights,  we  arranged  for  two  types  -- 
i)  those  including,  and  ii)  those  not  including  the  secondary  electron 
components  of  the  scattering  integral.  The  former  are  used  to  compute  the 
full  slowing-down  spectrum;  the  latter,  together  with  the  difference  between 
the  two  types,  are  necessary  for  computation  of  the  different  generations  of 
secondaries.  To  economize  memory  storage  space,  these  two  triangular  arrays 
of  integration  weights  are  both  stored  in  the  same  nearly  square  array 
(called  WT). 
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Table  II.l  gives  a summary  description  of  the  computer  program  which  we 
developed,  and  which  we  call  ELSPEL,  for  E^nergy  l^oss  ^ectrum  of  Electrons. 
Discussion  of  each  of  the  eighteen  routines  is  presented  in  Appendix  II. 

IV . Discussion  of  Results 

A.  Energy  Conservation 

We  begin  our  discussion  of  results  with  a table  showing  typical  energy 
and  number  conservation  data,  i.e.,  table  2.  A value  of  unity  is  the  desired 
result,  and  departures  from  unity  indicate  fractional  differences  from  source 
energies  T^.  Such  departures  are  seen  to  approach  a few  percent,  although 
for  most  energies  the  agreement  is  ~ 1%.  Since  the  integrals  required  for 
these  checks  feature  a large  contribution  from  the  region  of  T'  near  T, 
the  departures  from  unity  of  table  1 are  partly  due  to  evaluation  of  the  energy 
conservation  integral,  in  addition  to  errors  of  y.  But  for  the  most  part,  we 
believe  that  the  energy  conservation  integral  is  a reliable  measure.  The 
results  shown  are  sufficient  to  demonstrate  a general  feasibility  of. the 
approach  which  we  have  developed.  Further  improvements  in  computational 
accuracy  can  no  doubt  improve  this  check.  One  should  note  that  the  rate  of 
change  of  the  energy  conservation  data  is  probably  a better  indication  of  the 
accuracy  of  y in  any  given  part  of  the  spectrum  than  are  the  values  shown. 

Figure  5 gives  the  fractional  contribution  to  the  energy  deposition  by 
primaries,  secondaries,  and  all  other  components  combined.  Here,  energy 
deposition  includes  energy  losses  of  all  types  except  the  kinetic  energies  of 
secondary  electrons  still  capable  of  ionizing,  i.e.  kinetic  energies  > 16  eV. 
There  is  evident  a rather  remarkable  persistence  of  the  importance  of  the 
primaries  and  first  generation  secondaries,  which  together  are  respcfnsible 
for  more  than  90%  up  to  source  energies  greater  than  300  keV.  Since  both 
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components  can  be  represented  rather  well  with  fairly  simple  functional  forms, 
one  can  hope  to  obtain  analytic  representations  in  this  manner  which  are 
usable  for  many  types  of  yield  calculations. 

B.  The  Lewis  effect 

Figure  6 demonstrates  the  Lewis  effect  and  its  mai-^  features  at  a single 
(high)  source  energy  (1  MeV).  The  abscissa  measures  the  energy  loss  in  eV. 

The  Lyman  and  Werner  bands  appear  as  peaks,  with  multiple  interactions  producing 
a similar  structure  which  decreases  in  prominence.  The  maintenance  of  the 
similarity  of  structure  for  multiple  interactions  is  very  likely  due  to  the 
dominance  of  the  Lyman  peak:  probably  all  but  one  of  the  interactions  combining 

to  give  prominent  secondary  and  tertiary  peaks  are  Lyman  interactions. 

While  not  evident,  the  peaks  of  fig.  6 are  in  reality  composed  of  many 
individual  lines  which  are  not  resolved.  It  is  possible  to  resolve  these 
lines  with  finer  spacing,  birt  we  have  not  seen  fit  to  attempt  it. 

The  rather  arbitrary  cut-off  of  low-energy  knock-on  interactions  leads  to 
a rather  small  discontinuity,  labled  "Miller,"  which  is  identified  on  the 
figure.  The  form  chosen  for  this  cut-off  is  evidently  not  very  important. 

Also  displayed  is  the  energy-loss  region  for  which  energy  conservation  is  used 
to  correct  the  solution.  The  smoothness  apparent  in  this  region  is  partly  an 
artifact  of  the  method  used,  and  partly  a consequence  of  increased  interval 
size,  both  of  which  act  to  wash  out  the  physical  irregularities. 

The  ordinate  has  been  chosen  so  that  at  large  energy  losses  the  curve  in 
fig.  6 should  approach  unity,  corresponding  to  the  increasing  validity  of  the 
reciprocal  stopping  power  approximation  for  y(T,T^).  This  result  is  apparent 
although  the  asymptotic  trend  develops  more  slowly  at  larger  source  energies 
than  we  would  have  expected. 
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C.  Properties  of  si  owing-down  spectra 


Figure  7 presents  spectra  for  four  source  energies.  The  ordinate  has 
been  chosen  so  that  the  curve  is  fairly  flat  in  the  region  of  dominance  by 
secondary  electrons.^  For  the  30  keV  case,  contributions  by  primaries  and 
f i rst-generation  secondaries  are  indicated  by  dashed  lines;  and  in  the  case  of 
the  primaries,  the  dashed  line  has  approximately  the  functional  form  T/T^. 

According  to  the  CSDA  approximation,  the  first  generation  secondary 
electrons  due  to  a monoenergetic  source  at  are  described  by 

y(l)s  = f dT  K^(T',T)/S  , (9) 

"2T 

where  is  the  kernel  function  of  eq  (8).  If  we  write  the  stopping  power 

in  the  form  S = k(T)  B(T),  where  B(T)  is  the  stopping  number,  then  the 
lead  term  of  K^/S  is  [T~^  + (T'  - T)~"' ]/B(T' ) . We  assume  >>  T and 
make  use  of  the  usual  slow  rate  of  variation  of  B(T')  in  evaluating  this 
factor  at  T'  = T^/2,  which  is  approximately  the  midpoint  (T^  + 2T)/2  of 
the  integration  interval.  Then, 


S(T) 


y*'’(T,T 


1 oT  ^ T . ^o‘^ 

1 - 2—  + — log 


0 


T. 


0 


T 


T < T^/2 


1/B(y2),  T « T^. 


(9') 


The  minima  of  the  weighted  spectra  shown  in  fig.  7 are  thus  expected  to 
have  roughly  the  value  l/B(T^/2),  due  to  1 st-generation  secondaries  which 


^Since  the  product  S • y is  a rough  measure  of  the  number  of  electrons, 
(primaries  plus  secondaries)  per  primary,  at  any  given  kinetic  energy  T, 
the  ordinate  variable  can  be  interpreted  as  a rough  measure  of  the  fraction 
of  the  input  energy  T^  which  these  electrons  represent. 


22 


dominate  in  this  region.  Below  ^ 1 keV  the  minimum  tends  to  disappear. 

But  the  tendency  of  the  spectra  to  converge  keeps  a value  ~ B”''(.5  keV) 
significant  to  still  lower  energies.  At  the  right  of  the  figure,  some  of 
these  values  have  been  indicated. 

A remarkable  feature  noted  previously  by  Berger  [34]  is  a strong  tendency 
of  the  curves  to  converge  at  low  energies,  so  that  the  .:.hape  of  the  slowing- 
down  spectrum  is  nearly  independent  of  in  this  region.  One  sees  in 

eq  (9‘)  a tendency  of  this  type  in  the  y^^\  since  B varies  slowly  for 
large  T^.  In  higher  generations  of  secondaries  the  upper  limit  of  the 
integral  on  the  right  of  eq  (8)  becomes  largely  irrelevant,  so  that  one  expects 
these  components  to  depend  on  only  through  B(T^/2)  at  most.  But  the 

convergence  seen  in  fig.  7 is  more  striking  than  this  and  must  reflect  further 
restrictions  imposed  by  energy  conservation.  Table  3 illustrates  this  strong 
trend  towards  T^-independence  of  the  low  energy  components.  Convergence  is 
evidently  according  to  a trend  somewhat  of  the  type 

y(T,T^)  ~ y(T,-)  + a^/log(T^/a^)  , (10) 

where  a^,  a-j  depend  on  T only,  because  equal  factors  of  increase  in  source 
energy  improve  convergence  by  comparable  i ncrements . 

D.  The  ionization  yield 

Because  degradation  spectra  have  extensive  applications  to  the  evaluation 
of  excitation  and  ionization  yields,  we  turn  to  the  discussion  of  the  ionization 
yield,  usually  given  in  terms  of  a quantity  W,  eV  per  ion  pair, 

= ty’  / ^i(T)  y(T,y)  , 

0 

where  k^.  is  probability  for  ionization  per  unit  path  length  by  electrons 
with  kinetic  energy  T.  At  the  bottom  of  fig.  7 a curve  of  K_j/S  has  been 
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included.  It  is  instructive  to  consider  the  different  contributions  to  the 
value  of  W"^  for  H^.  The  T/T^  trend  of  fig.  7 for  the  primaries,  which 
follows  from  CSDA,  has  the  consequence  that  the  contribution  by  primaries  can 
be  approximated  by  evaluating  (k^./S)  at  the  midpoint  energy  T^/2,  e.g., 

.018  for  Tq  = 30  keV.  Correspondingly , the  fact  that  both  (k^./S)  and 
(T/T^)Sy  are  comparatively  flat  means  that  a rough  evaluation  of  the  contri- 
bution due  to  secondaries  can  be  obtained  with  similar  ease.  One  notes  that 
large  contributions,  say  in  the  case  of  30  keV,  occur  between  ~ 20  eV,  and 
~ 2500  eV.  The  midpoint  of  this  range  (in  the  logarithm)  is  (20  x 2500)^  = 

224  eV.  Evaluation  of  both  K^/S  and  (T/T^)Sy  curves  at  this  energy,  with 
the  product  multiplied  by  £n(2500/20)  = 4.83  gives  a value  .0095.  We  thus 
estimate  that  W 1/(.018  + .0095)  36.4  eV.  This  agrees  rather  well  with 

more  accurate  calculations. 

One  can  see  from  these  curves  that  the  high  energy  gradual  decrease  of 
k^/S  will  decrease  the  contribution  to  W~^  due  to  primaries,  as  the  source 
energy  T^  is  raised.  On  the  other  hand,  secondary  electron  contributions 
will  increase.  Unfortunately,  the  resulting  rough  cancellation  of  competing 
effects  provides  no  simple  explanation  for  the  impressive  T^-independence 
found  in  experimental  investigations  of  W.^  But  this  approach  could  probably 
be  developed  to  give  rough  estimates  fairly  readily. 

Pursuing  the  question  of  W further,  we  present  values  due  to  other 
investigators  in  table  4,  and  for  different  source  energies  in  table  5.  A 
minimum  occurs  in  the  vicinity  of  2 keV,  with  a very  slow  rise  with  increasing 
source  energy.  Comparisons  with  results  from  experiments  and  other  recent 

6 

See  refs  (7)  and  (8)  for  an  interesting  approach  to  this  problem  which  uses 

scaling  of  the  si  owing-down  spectra. 
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computations  are  very  interesting,  and  are  displayed  in  figs.  8 and  9.  In 
fig.  9 one  sees  that  Garvey,  Porter,  and  Green  (GPG)  obtain  a much  larger  rate 
of  increase  in  W with  source  energy  above  ~ 3 keV.  On  the  other  hand,  in 
fig.  8 the  data  of  GPG  show  a more  rapid  increase  of  W as  T is 
decreased  from  3 keV  to  below  200  eV.  We  have  not  attempted  to  analyze 
these  differences  in  detail;  but  analogous  differences  in  the  case  of  H2O 
can  apparently  be  traced  to  differences  of  cross  section.  One  need  not 
assume  that  computational  procedures  are  the  cause  of  such  contrasting  trends. 

The  case  is  different  in  the  comparison  with  Gerhart  shown  in  fig.  8. 

While  the  cross  sections  used  in  our  calculations  are  not  identical  with  those 
of  Gerhart,  the  differences  would  not  seem  to  be  sufficiently  great  to  give 
rise  to  the  ~ 10%  differences  in  W which  occur. 

In  fig.  8,  the  dashed  line  gives  values  of  W inferred  from  the  following 
formula,  which  depends  on  the  ratio  of  total  to  ionization  cross  sections: 


W 


- ^tot  |~^o 
° ^Toin^^ 


(11) 


As  discussed  by  Jones  [51],  such  an  expression  is  valid  when  one  ionization  is 
energetically  possible,  but  two  interactions  of  which  at  least  one  is  an 
ionization  are  not  possible.  This  condition  applies  below  about  26  eV  for 
the  case  of  H^.  The  bump  shown  in  the  W curve  is  therefore  a consequence 
of  cross  section  ratios  and  does  not  follow  from  multiple  processes. 
Experimental  values  due  to  Combecher  [57]  do  not  show  this  reversal  clearly, 
although  they  exhibit  a shoulder  at  correspond! ng  energies.  Such  reversals 
have  been  reported  in  references  [23,58,22],  and  also  by  D.  Douthat  [54].  We 
do  not  know  whether  this  is  a real  phenomenon  or  not.  It  is  strongly 

C 

dependent  on  cross  sections  which  are  not  well  known  as  yet. 
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Equation  (11)  furnishes  a good  monitor  for  the  accuracy  of  the  slowing- 
down  computations,  since  corresponding  values  can  be  obtained  by  application 
of  the  computer  programs,  which  were  designed  for  much  greater  source  energies 
and  are  therefore  put  to  a hard  test  at  low  energies.  Agreement  shown  in 
fig.  8 thus  constitutes  a good  check. 

Figure  10  gives  a different  plot  of  some  of  the  data  of  fig.  8.  Below 
~ 27  eV,  the  solid  curve  represents  the  cross  section  ratio  of  eq  (11). 

Below  20  eV  Combecher's  data  differs  by  a factor  ~ 1.2  in  magnitude  from  the 
solid  curve;  and  by  inference  from  eq  (11),  this  difference  measures  differences 
of  cross  sections.  Since  our  calculations  thus  appear  to  permit  too  much 
ionization  at  low  energies  our  W values  may  be  somewhat  low  at  high  energies. 
But  this  same  argument  means  that  Douthat's  two  versions  of  cross  sections 
below  100  eV  should  give  rather  different  W values  at  high  energies.^  In 
our  experience,  high  energy  W values  are  rather  sensitive  to  low  energy 
cross  sections,  even  those  below  25  eV. 

It  may  be  that  the  rather  reasonable  values  for  W obtained  here  occur 
because  Gerhart's  cross  sections  for  ionization  are  too  large  at  low  energies, 
as  evidenced  by  Combecher,  but  are  correspondingly  too  small  at  high  energies, 
as  much  other  data  suggests. 

Figure  11  compares  our  slowing-down  spectra  for  = 10  keV  with 
corresponding  spectra  due  to  Douthat  [1]  and  GPG  [2],  as  read  from  the 
published  curves.  In  a general  way,  our  values  lie  between  those  of  the  other 
authors,  with  the  largest  differences  occurring  at  low  energies,  where  factors 
of  2 can  be  seen.  For  the  region  above  ~ 300  eV,  which  is  dominated  by 
y - S \ the  agreement  is  generally  good,  as  can  be  expected. 

^Douthat's  lower  values,  € in  fig.  7,  correspond  to  a cross  section  set 
rather  similar  to  ours.  His  higher  values,  +,  correspond  to  about  a 28% 
reduction  in  excitation  cross  section  with  16%  increase  in  ionization 
cross  sections. 
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100  keV  due  to 


Table  6 gives  a few  comparisons  with  data  for  = 

Okazaki  and  Sato  [15].  No  agreement  is  evident  even  at  high  energies  where 
the  y values  depend  essentially  on  the  stopping  power  alone. 
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Appendix  I.  Cross  Section  Data 


For  purposes  of  constructing  a computer  program  for  the  integration,  we 
have  used  a well-known  set  of  cross  sections,  viz  the  cross  section  set 
developed  by  Donald  E.  Gerhart  for  his  PhD  thesis  [3],  even  though  there  are 
some  unsati sfactory  aspects  to  them.  In  this  'appendix  we  describe  necessary 
modifications  and  adaptations. 

Gerhart  uses  Platzman's  approach,  which  depends  on  the  Wil 1 iams-Wei szacker 
approximation  in  which  the  passage  of  a fast  charged  particle  is  replaced  by  a 
nearly  equivalent  spectrum  of  virtual  photons,  which  are  subject  to  photo- 
electric absorption  [84].  Collisions  are  classified  as  either  "glancing''  or 
"knock-on."  Glancing  collisions,  with  small  momentum  transfer,  are  assumed  to 
follow  expressions  derived  for  zero  momentum  transfer.  Knock-on  collisions 
are  assumed  to  follow  cross  section  forms  based  on  collisions  between  free 
el ectrons. 

A.  Glancing  Collisions  (Small  momentum  transfer,  large  impact  parameter). 

First,  in  regards  to  excitations  due  to  glancing  collisions  we  follow 
Gerhart  in  the  use  of  the  Bethe  excitation  cross  section  formula  as  follows: 


Aira^R  f 

0 n „ , 4T  N 

T E^/R  R ^set^  ’ 


(I.l) 


where  R = 13.6  eV  and  f^  is  the  oscillator  strength  for  the  n'th 
transition  of  a set,  which  can  be  any  group  of  transitions  treated  together. 

The  constant  c^^^  in  the  logarithm,  which  in  principle  differs  from  transition 
to  transition,  is  assigned  an  average  value  appropriate  to  a set  of  related 
lines  such  as  those  constituting  the  Lyman  or  Werner  bands  of  ^ 

Ionizations  due  to  glancing  collisions  are  described  by  the  analogous 
cross  section  form 
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a . 


(1.2) 


where  is  an  "ionization  efficiency"  function,  and  b is  essentially 

the  reciprocal  of  an  impact  parameter  in  a^  units,  or  a cut-off  for  the 
square  of  the  momentum  transfer  in  glancing  collisions.  The  variable  E 
is  the  transition  energy,  and  df/dE  is  the  differential  oscillator  strength 
distribution  (OSD). 

The  OSD  developed  by  Gerhart  is  based  largely  on  discrete  excitation 
data  due  to  Allison  and  Dalgarno  [59]  and  experimental  data  due  to  Cook  and 
Metzger  [60],  both  modified  to  satisfy  information  from  sum  rules.  The 
important  and  complex  far-ul traviol et  region  from  14.4  to  16.5  eV  is  described 
by  replacing  the  discrete  Rydberg  states  with  a carefully  estimated  continuous 
extension  of  the  continuous  part  of  the  OSD  at  higher  energies. 

The  f^  data  of  Table  I.l  are  from  Gerhart,  and  represent  his  modifica- 
tion of  the  data  due  to  Allison  and  Dalgarno.  The  corresponding  E^  data 
were  not  given  nor  referred  to  in  [3];  they  are  due  to  Namioka  [29]. 

The  ionization  efficiency  function  n.j  used  by  Gerhart  is  reproduced 
from  his  paper  in  fig  I.l.  A dash-dot  line,  corresponding  to  an  analytic 
approximating  function 


has  been  added;  we  used  this  curve  for  convenience,  despite  the  fact  that  it 
can  clearly  be  improved. 

The  mathematical  form  of  eq  (I.l)  is  a discrete  analog  of  eq  (1.2).  In 
both  cases,  the  expression  can  be  written  as  the  sum  of  two  products,  each 
having  a function  f(T)  multiplied  by  a function  g(E)  (or  g^)  of  the 
transition  energy.  We  proceed  to  a unified  treatment  of  excitations  and 


(1.3) 
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ionizations  in  which  cumulative  integrals  (or  sums),  extend  from  zero  transition 
energy  over  excitations  and  then  ionizations,  up  to  a maximum  energy  of  about 
220  eV.  Three  such  cumulative  integrals  (f^,  ^ multiplied  by  factors 
or  e\  i = 0,1,2),  are  tabulated  at  intervals  of  about  0.4  eV  for  small 
transition  energies,  with  the  interval  size  increasing  geometrically  with  E. 

A second  set  of  integrals  includes  a factor  £n(E^)  or  £n(E).  The  resulting 
6 tabulations  are  given  in  Table  1.2.  The  first  three  columns  {£)  are  energy 

j 

moments  of  f^  and  while  the  latter  three  columns  (^ ) are  energy 
moments  of  ^^n(E),  E^  and  E in  eV.  Factors  (1  - n^.) 

and  n.j  are  included  in  calculations  for  excitations  and  ionizations, 
respectively.  These  tables  form  the  numerical  basis  for  all  computations 
involving  glancing  collisions. 

The  use  of  three  integrals,  yielding  three  moments  of  the  transition 
energy,  has  the  following  significance;  The  zeroth  momenth  is  essential  to 
the  evaluation  of  the  total  cross  section;  the  first  moment  is  required-  for 
stopping  power  calculations;  and  by  use  of  three  moments  one  may  perform  Gauss 
quadrature  calculations  of  quadratic  accuracy.^ 

The  constant  b was  assigned  a value  2.3  by  Gerhart.  As  will  be 
discussed  in  section  C below,  we  found  it  convenient  to  make  b a function 
b(T)  for  electron  kinetic  energies  T below  707  eV. 

We  note  that  eqs  (I.l)  and  (1.2)  are  non-relativistic  and  require 
modification  to  achieve  the  relativistic  forms.  This  is  essential  because  we 
wish  to  do  calculations  at  relativistic  as  well  as  non-relativistic  e'nergies. 
Hence  we  modified  all  glancing  cross  sections  by  replacing  the  logarithmic 
factor  ln(4T/R)  as  follows: 


^The  second  moment  is  important  in  straggling  theory. 
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->  £n 


(T 


+ 2) 


* * 


T (T  +2) 
(T  +1)^ 


(1.4) 


•k  2 

Where  T = (T/mc  ) and  the  term  on  the  right  is  the  square  of  the  velocity, 

2 2 

3 = (v/c)  . Likewise,  the  factor  (R/T)  appearing  at  the  left  of  cross 

2 2 

section  forms  given  in  eqs  (I.l)  and  (1.2)  was  replaced  by  2(R/mc”)3  . 

B.  Knock-on  Collisions  (large  momentum  transfer,  small  impact  parameter). 
Gerhart  also  made  use  of  non-relativistic  forms  for  the  knock-on  cross 
section.  These  he  multiplied  by  a factor  c{)(E)  and  by  another  factor  C(T). 

The  latter  factor,  which  was  evidently  introduced  to  improve  agreement  with 
total  cross  sections  below  ~ 1000  eV,  makes  it  impossible  to  derive  the  Bethe 
stopping  power  formula  without  some  change  in  Gerhart's  formulae.  Hence  we 
simply  eliminate  C(T)  by  equating  it  to  unity.  This  leaves  us  with  Gerhart's 
problem  of  obtaining  agreement  with  data  on  total  cross  sections  below  1 keV, 
and  this  is  discussed  in  the  next  section. 

Gerhart's  function  (f)(E)  was  introduced  as  a smooth  cut-off  to  the 
spectrum  of  transition  energies  resulting  from  the  knock-on  process.  (1)(E) 
is  shown  in  fig.  1.2,  together  with  a simple  step-function  approximation  which 
assumes  zero  cross  section  for  E < 17.5  eV,  which  we  decided  to  use.  In 
our  opinion,  the  slowing-down  spectra  are  not  sensitive  to  such  a change  of 
shape;  and  the  sharp  cut-off  simplified  both  numerical  and  analytic  procedures. 

Consistent  with  the  point  of  view  already  expressed,  we  use  the  fully- 
relativistic  M(|)ller  cross  section  to  represent  the  knock-on  process,  modified 
in  a simple  way  to  take  account  of  a mean  energy  I « 1 6 eV  for  ionization: 
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+ 


(1.5) 


+ 


u(E-17.5)  , 


where  the  unit  function,  u(x)  = 


1 , X > 1 , 

0,  X < 1 , 


expresses  the  cut-off  shown  in 


fig.  1.2  for  low  energy  excitations  by  this  process. 


C . Low  energy  ionization  cross  sections. 

_2 

Gerhart  chose  to  use  a simple  E cross  section  for  all  processes  at 


low  energies,  i.e.,  T < 300  eV.  This  function  was  normalized  to  give  a value 
for  the  total  ionization  cross  section  in  agreement  with  experimental  data  due 
to  Rapp  and  Englander-Golden  (REG)  [33].  We  have  chosen  a different  approach 
to  the  problem  of  the  selection  of  differential  cross  sections  at  low  energies: 
we  continue  to  use  eq  (1.2),  but  assign  to  b (for  T < 707  eV)  an  energy 
dependence  b(T)  adjusted  to  give  total  cross  sections  in  agreement  with 
REG.  The  theoretical  basis  for  the  T-independence  of  b breaks  down  at  low 
energies.  In  fact,  since  there  is  no  theoretical  basis  at  low  energies  for 
the  asymptotic  form  which  we  use,  this  procedure  should  be  considered  only  as 
an  ad  hoc  extrapolation  guided  by  an  experimental ly  fixed  integral. 
Incidentally,  Miller  used  b = b(E)  [58],  but  there  is  some  theoretical 
justification  for  this  [62].  Even  so,  as  the  glancing  component  decreases  in 
importance  at  low  T,  it  is  not  unreasonable  to  assign  this  to  a decreasing 
range  of  contributing  momentum  transfers. 

Use  of  values  for  b(T)  as  large  as  b = 3 at  T = 1 keV  could 
accomplish  essentially  the  same  thing  as  Gerhart's  C(T)  although  we  have  not 
incorporated  this  adjustment;  while  for  T < 500  eV  use  of  b(T)  results  in 
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more  gradual  change  of  spectrum  than  that  due  to  Gerhart's  prescri ption . 

Table  1.3  gives  £n[b(T)]  values  which  we  have  used.  Table  1.4  shows  the 
comparison  with  REG  thus  attained. 

Our  approach  is  no  less  arbitrary  than  that  adopted  by  Gerhart  and  in 
modified  form  by  Douthat.  Our  value  for  b above  707  eV  is  2.34.  This 
sacrifices  close  agreement  with  REG  in  the  region  707  eV  < T < 1000  eV,  but 
preserves  the  asymptotic  agreement  of  stopping  powers  at  high  energy  with  the 
Bethe  stopping  powers  based  on  a value  for  the  mean  excitation  potential 
consistent  with  Gerhart's  OSD. 

D.  Low  energy  excitation  cross  sections 

Figures  9 and  10  of  Gerhart's  paper  [3]  present  his  "adopted"  electron 
impact  cross  sections.  Since  our  purpose  here  is  to  perform  electron  slowing- 
down  calculations  with  minimal  changes  in  Gerhart's  cross  section  set,  we 
have  reduced  these  curves  to  tabular  form  for  our  use,  with  no  significant 
modifications.  But  one  should  note  that  accurate  values  in  the  important 
energy  region  13  eV  < T < 25  eV  cannot  be  obtained  in  this  way,  and  detailed 
features  of  the  spectra  in  this  region  can  be  wholly  spurious. 

These  total  cross  section  data  must  be  in  some  manner  complemented  with 
information  about  associated  transition  energies.  In  the  case  of  the  Lyman 
and  Werner  discrete  excitations  one  has  values  for  individual  transition 
energies  [29];  but  for  higher  energy  excitations  a range  of  energies  from 
approximately  14.4  eV  to  17  eV  is  involved,  with  the  use  of  a continuous  OSD 
function  in  place  of  the  actual  complex  line  structure  as  already  noted.  Our 
first  approach  was  to  choose  for  each  process  separately  a functional  form  of 
the  type 

[1  - n,(E)]  I ^^n(4lb(T)  . (1.6) 
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But  since  the  si  owing-down  calculation  depends  for  the  most  part  only  on  three 
moments  (integrals)  over  the  transition  energy,  we  carried  use  of  this  form  a 
step  further  and  applied  it  to  the  combination  of  al 1 excitation  cross 
sections.  This  leads  to  slightly  different  values  b.(T),  j = 1,2,3,  for  the 

vj 

three  different  moments,  as  required  to  agree  with  calculations  having  different 
b(T)  values  for  different  processes  as  in  eq  (1.6).  However,  it  turned  out 
that  these  b-  differed  in  a regular  manner,  but  only  by  ~ 1%.  Hence  we 

<J 

eventually  used  only  b-j(T),  which  was  assigned  values  very  slightly  larger 
than  those  calculated  from  total  cross  section  data.  Tabulated  data  for 
resulting  b-j  function  are  given  in  Table  1.3;  and  a comparison  of  resulting 
total  cross  section  values  to  those  of  Gerhart  is  given  in  Table  1.4. 

The  full  differential  spectrum  is  important  when  single  interactions  are 
significant.  We  calculated  these  spectra  by  use  of  eqs  (1.2)  and  (1.6),  using 
values  for  £n(b)  and  £n(b-j)  already  described.  Figures  1 and  2 give  some 
results  of  these  calculations. 

Vibrational  effects  become  important  below  about  20  eV;  and  they 
constitute  the  main  part  of  the  cross  section  below  the  lowest  electronic 
transition  energy,  here  taken  to  be  11.213  eV.  Though  ignored  by  Gerhart, 
we  have  included  them  on  the  basis  of  data  by  Ehrhardt  et  al . [63]  and 
Trajmar  et  al.  [64].  We  use  only  two  excitation  energies,  0.54  eV  and 
1.08  eV.  For  electrons  with  kinetic  energy  greater  than  20  eV,  we  use  the 
following  simple  expressions  for  the  (non-relativistic)  differential  cross 
sections : 

k^/K  [-  .00457  + .00588  £n(4T/R)]  5(E  - .54)  , (1.7) 

k2/K  ^ [-  .000309  + .000287  £n(4T/R)]  6(E  - 1.08)  , (1.8) 
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where  k = 4iTa^R/T.  Below  20  eV,  we  tabulated  the  total  cross  section  data 
and  perform  interpolations  as  required.  Tabular  data  are  recorded  in 
Table  1.5;  these  are  a somewhat  smoothed  version  of  the  data  in  Ehrhardt  et  al 
and  Trajmar  et  al . 

The  experimental  data  do  not  e: fend  beyond  80  eV;  hence  the  use  of 
eqs  (1.7)  and  (1.8)  for  much  higher  energies  is  not  necessarily  accurate.  But 
because  the  contribution  is  a small  one,  this  extrapolation  was  used  (with  the 
relativistic  modifications.)  Equations  (1.7)  and  (1.8)  are  of  the  type  which 
can  be  readily  incorporated  into  the  subroutine  OSCUM  which  evaluates 
differential  and  integral  total  oscillator  strengths.  But  we  prefer  to  treat 
vibrations  by  use  of  separate  computations  which  could  be  used  or  not  as 
desi red. 

E.  Symmetry  of  glancing  collision  cross  sections 

One  problem  which  must  be  solved  to  preserve  energy  conservation  at  low 
energies  is  the  following:  In  the  case  of  glancing  collisions,  for  which  the 

cross  sections  are  not  written  symmetrical ly , secondary  electrons  with  kinetic 
energy  (E  - I)  are  generated  with  the  same  probability  as  that  corresponding 
to  energy  losses  E.  A type  of  symmetrizing  is  necessary  in  order  to 
guarantee  this  feature.  We  chose  to  accomplish  this  by  terminating  the  cross 
sections  at  an  energy  loss  for  which  secondary  and  primary  electrons  had  equal 
kinetic  energy,  and  then  repeating  the  cross  section  structure  below  this 
energy  so  that  it  has  identical  values  for  kinetic  energies  equally  above  and 
below  that  value  at  which  primary  and  secondary  kinetic  energies  coincide. 

This  seems  a little  easier  than  the  alternative  approach  of  adding  a 
symmetrizing  term.  Where  the  difference  between  approaches  is  significant,  we 
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feel  that  the  effect  on  the  spectrum  is  small,  so  that  any  improvement 
would  call  for  evaluation  of  more  accurate  and  properly  symmetrized  cross 
sections,  i.e.,  a major  effort. 

F.  The  mean  excitation  potential 

A mean  excitation  potential  <I)g  can  be  determined  for  the  oscillator 
strength  distribution  for  use  in  Bethe's  stopping  power  formula  [65]: 


<I)n  = 2 R exp  SZ' 


r 


df  £n  (E/2R) 


] 


(1.9) 


Analogously,  cross  sections  which  we  have  adopted  lead  to  a high  energy 
stopping  power  of  Bethe's  form;  and  if  we  determine  an  equivalent  parameter 
<I)p  (P  for  Platzman)  we  find  the  following  expression: 


<I)p  = 2R  exp-jz‘ 


j df  £n 
*-noniz. 


/bR 


^ n ^ set' 
sets  n 


+ 2 (4R 


(I. 10) 


where  is  a cut-off  for  low  energy  knock-on  collisions,  which  was  assigned 

the  value  17.5  eV  to  correspond  to  Gerhart's  cut-off  function  j)(E). 

Table  1.6  gives  values  for  <I)g  and  <I)p,  in  comparison  with  values 
by  other  authors.  Consistency  between  <I)g  and  <I)p  can  be  improved,  but 
it  is  interesting  that  (I)g  agrees  rather  well  with  other  theoretical 
estimates,  while  <I)p  is  consistent  with  experimental  values.  Partly  for 
this  reason  we  have  not  enforced  agreement  between  <I)g  and  <I)p. 
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APPENDIX  II.  Description  of  the  Subprograms  constituting  ELSPEL 


The  descriptions  which  follow  are  in  the  order  shown  in  table  II. 1, 
proceeding  from  differential  to  integral  cross  sections,  to  the  subprograms 
guiding  the  solution. 


1.  OSCUM(E,ISET,ID) 

The  output  from  this  subroutine  consists  either  of  cumulative  integrals 
over  the  oscillator  strength  (OS)  distribution  (ID  =0),  or  the  differential 
distribution  itself  (ID  =1).  The  six  cumulative  integrals  (for  ionization) 
are  as  follows:  (R  = 13.6  eV) 


where  E is  in  eV  in  the  factor  £n(E). 

For  excitations  the  corresponding  quantities  are,  (E^  e AE)  and  m = 0,1,2, 


A.  Elementary  Cross  Section  Subroutines 


(II. 1) 


0 


m-1  ^set 
n 


(II. 2) 


E E(E„/R) 

sets  n 


m-1  ^set 
n 
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where  the  distinct  "sets"  include  Lyman,  Werner,  excitations  with  principle 
quantum  number  n = 3,  and  the  combination  of  all  other  excitations  involving 
transition  energies  between  14.2  eV  and  16  eV.  But  we  reduce  these  four 
types  of  interaction  to  three  in  practice  because  the  OS  distribution,  which 
is  taken  to  be  effectively  continuous  in  this  energy  region,  is  treated  in 
Gerhart's  paper  as  identically  the  same  in  shape  for  all  excitations  which  are 
neither  Lyman  nor  Werner  [3]. 

Cumulative  integrals  over  excitation  oscillator  strength  appear  in  an 
output  array  ASRE(I),  where  I = 1,2,3  correspond  to  m = 0,1,2  of  the 
integrals  of  eqs  (II. 2),  and  I = 4,5,6  relate  to  the  corresponding 
integrals  of  eqs  (II. 2).  Similarly,  for  ionizations  the  output  is  in  ASRI(I), 
where  index  values  I = 1,2,3  correspond  to  the  first  integrals  of  eqs  (II.l); 
and  index  values  I = 4,5,6  correspond  to  the  integrals  of  eqs  (II.l) 
containing  the  logarithmic  factor. 

For  ID  = 1 and  E ^ 17  eV,  the  value  of  df/dE  for  ionizations 
appears  in  ASRI(l),  while  £n(E)df/dE  appears  in  ASRI(2).  For  E < 17  eV, 
equivalent  quantities  are  obtained  by  differencing  the  cumulative  tabulations; 
these  appear  in  ASRE(l),  ASRE(2),  ASRI(l),  and  ASRI(2),  the  index  2 
referring  to  the  quantity  with  the  logarithmic  factor.  (These  equivalent 
quantities  are  and  divided  by  the  width  AE  for  any  interval.) 

Note  that  for  the  ID  = 1 case,  the  output  has  units  (eV)~"'. 

Input  tabulations  for  OSCUM  are  read  from  cards,  when  ISET  f 1976. 

One  set  of  input  data  (EON,  for  ionizations)  consists  of  cumulative  'integrals 
tabulated  at  intervals  from  15.176  eV  to  17  eV.  The  other  set  (EXC,  for 
excitations)  is  a similar  tabulation  extending  from  11.213  eV  to  17  eV. 

The  two  tabulations  are  listed  in  table  1.2.  Correspond! ng  energies  are  read 
into  EOS. 
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The  method  for  computing  ASRE  and  ASRI  values  is  interpolation  on  the 
prepared  tables.  For  our  calculations,  based  on  Gerhart's  set  of  cross  sections, 
the  following  forms  were  used  to  extend  the  tables  to  222.7  eV,  where  E 
values  are  in  eV: 

transition  energy  (list  name:  EOS): 


En  = E^^_^  (17/11. 213)^/^\  n 


. . . , 80 


(II. 3) 


excitations  and  ionizations,  14.2  eV  to  17  eV  (list  name:  EXC): 


(II. 4) 


ionization  efficiency: 


(II. 5) 


ionizations  (list  name:  EON),  17  - 35  eV: 


^ ^ 275 


(II. 6) 


dE  ^2.77 


for  35  eV  to  222.7  eV: 


df  1624 
dE  ^3.27 


(II. 7) 


for  E > 222.7  eV 


(II. 8) 
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Interpolation  requires  both  positioning  on  the  table,  and  an  interpolation 
formula.  Positioning  was  accomplished  as  follows  by  determining  the  equivalent 
index  by  the  formula  below,  where  N is  the  next  integer  below  EX: 


PY  - loq(E/11.213) 

(1/11 )log(17/ll .213) 


(II. 9) 


N 4 EX 


(II. 9') 


The  interpolation  weights  could  then  be  calculated  from  (EX  - N),  with 
integer  separation  between  tabulated  abscissa  values.  A simple  3-point  rule 
was  used  for  this  purpose. 

The  df/dE  tabular  values  were  obtained  as  previously  mentioned  by 
differencing  as  follows  in  the  region  below  17  eV: 

E0N(I,J+1)  - E0N(I,J) 

E0S(J+1)  - EOS(J)  ’ 


EXC(I,J+1)  - EXC(I,J) 
E0S(J+1)  - EOS(J) 


for  1=1,4  and  1 < J < 12. 


(II. 10) 


No  attempt  was  made  to  apply  interpolation  for  E < 17  eV  when  df/dE  was 
required;  a histogram  structure  was  assumed.  Interpolations  were,  however, 
used  for  ionization  above  17  eV. 


2.  VIB(T,ARGL,X1 ,X2) 

Output  for  this  subroutine  are  values  for  the  total  vibrational  cross 
sections  at  kinetic  energy  T,  for  1st  and  2nd  excitation  levels;  these 
values  are  designated  XI  and  X2,  respectively,  for  energies  greater  than 
T = 20  eV,  the  values  XI,  X2  are  given  by  formulae  as  follows: 
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XI 


.00457  + .00588  * ARGL 


X2  = - .000309  + .000287  * ARGL,  (1 1. 11) 

where  the  constants  .00588  • (.54/R)  and  .000287  • (1.08/R)  are  oscillator 
strengths  for  the  two  transitions. 

In  adaption  to  the  input  parameter 

ARGL  = £n[2(T/R)(2  + T/mc^)]  - 6^(T)  , (11.12) 

an  i nput  table  VB  for  each  level  is  read  from  cards,  as  is  a corresponding 
tabulation  of  the  quantity  iln(4T/R),  which  is  stored  in  the  list  TVB. 

For  energies  less  than  20  eV,  the  method  used  is  an  interpolation,  by 
means  of  3-point  weights.  Positioning  on  the  TVB  table  for  a calculation 
at  T = is  accomplished  by  a search  for  the  lowest  tabular  value  exceeded 
by  U = X,n(4T,/R) . 

a 

Data  on  which  the  linear  forms  (for  T > 20  eV)  are  based  do  not 
extend  above  80  eV;  hence  the  constants  used  are  probably  not  very 
accurate  at  higher  energies.  But  the  cross  sections  are  small  at  high 
energies,  so  that  this  is  not  important. 

3.  BM0D(T,B1) 

The  output  value  of  this  function  routine  is  B1  and 

BMOD  = X^n(bjb)  , (11.13) 

where  b = b(T)  is  the  momentum  parameter  which  separates  glancing  from 
knock-on  collisions.  For  the  tabulation  we  have  taken  b = b^  for 
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T > (1000)/2  eV,  and  have  used  ^ 3.  For  lower  electron  kinetic  energies, 
the  value  of  b(T)  has  been  chosen  to  render  a total  cross  section  for 
ionization  in  agreement  with  experiment  (Rapp,  Englander-Golden  [29]).  But 
in  practice  we  modify  this  to  £n(b^/b)  = .25  for  T > 707  eV,  and  b^ 


correspondi ngly  equals 

2.34. 

The  parameter  B1 

is  defined  similarly, 

B1  = £n[b-|  (°°)/b-|  ] (11.14) 

where  b-j  (°°)  = 1,  a value  which  applies  for  T > 350  eV,  as  shown  in  Table  1.3 
This  parameter  corrects  the  total  cross  sections  to  agree  with  low  energy 
values,  as  discussed  in  section  3(D). 

Input  consists  of  a tabulation  of  values  for  £n(b^/b)  and  i^n[b.j  (“)/b-|  ] , 
which  are  read  from  cards  into  the  lists  BC0R(I,1)  and  BC0R(I,2),  respec- 
tively. Data  are  given  in  Table  1.3. 

Method  of  evaluation  is  interpolation  on  the  table.  Between  15.62  eV 
and  707  eV,  2-point  trapezoidal  rule  interpolation  is  used.  Below  15.62  eV, 
a 3-point  rule  is  used  instead.  Positioning  on  the  table  is  accomplished 


as  in  OSCUM,  viz.  for 

N equaling  the  integer  just  less  than  X, 

V _ loq(1000/T)  , , 
log/? 

(II. 15) 

N < X 


An  extra  point  has  been 
accuracy  near  threshold. 

inserted  at  low  energies  in  each  list  to  give  greater 
, ^ 
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4.  DTU(T,B,TB) 

Output  of  this  function  routine  are  two  numbers,  DTU  and  TB.  Both 


relate  to  maximum  transition  energy  for  glancing  ionizations.  For  a given 
value  of  T,  and  corresponding  cut-off  b(T),  DTU  is  the  maximum  allowed 
energy  for  glancing  ionization  transitions  from  higher  kinetic  energies  T. 
Correspondingly,  TB  is  the  least  kinetic  energy  for  electrons  capable  of 
giving  rise  to  the  transition  energy  (T  + I),  where  T is  the  kinetic 
energy  of  secondary  electrons  and  I = 16  eV  is  a mean  ionization  energy  for 
molecular  hydrogen.  Thus  DTU  relates  to  the  energy  loss  of  primaries, 
and  TB  to  the  energy  of  secondaries;  both  express  glancing  collision  limits. 

In  addition  to  T and  B = £n(b^/b(T)),  this  subroutine  requires  BMOD , 
because  both  DTU  and  TB  are  evaluated  by  reverse  interpolation  to  solve 
the  transcendental  equations 

* * / T 

2T  Rb(T  ) 2 + ^ 

\ mc"^ 

(11.16) 

★ * / T* 

2T  Rb(T  ) 2 + -!-2 

\ mc'^ 

The  reverse  interpolation  is  based  on  a list  of  energies 

RAD(I,1)  = 9(2^/^)^'^  eV,  (11.17) 

and  a corresponding  tabulation  RAD(I,2)  of  values  for  the  square  root  of  the 

left  side  of  the  transcendental  equations  above,  when  values  of  the  energy 

★ 

list  RAD(I,1)  are  equated  to  T . 


exp[-B^(T*)]  = (T  + I)^  , T*  = TB 


j exp[-B^(T*)]  = (DTU)^  T*  = T + DTU 
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The  positioning  is  accomplished  by  a search  to  locate  both  (T  - T) 
and  (T  + I)  in  the  list  RAD(I,2),  treated  as  abscissa  list.  A 3-point 
interpolation  is  then  applied  to  RAD(I,1),  treated  as  the  ordinate  list,  to 

•k  -k 

find  the  value  of  T corresponding  more  exactly  to  (T  - T)  and  (T  + I) 

•k 

than  the  neighboring  tabulated  values  of  RAD(I,2).  This  T then  gives  the 

•k 

output,  directly  in  the  case  of  TB,  and  through  use  with  (T  - T)  in  the 
case  of  DTU. 

This  rather  lengthy  positioning  operation  would  probably  take  too  much 
time  if  DTU  and  TB  values  were  required  throughout,  rather  than  for  the 
list  of  solution  energies  only,  i.e.,  a linear  set  of  no  more  than  100  to  200 
values. 


B.  Combined  Cross  Section  Functions 

5.  TCX(T,DT,IVB) 

Output  of  this  function  routine  is  a quantity  which  can  be  considered  as 

2 

me  ^ 2 2 

either  (for  non-relativistic  energies)  total  cross  section  times  (T/Arra^R'^) 
or  total  interaction  probability  divided  by  XKAP  (as  defined  below),  for  all 
processes  at  kinetic  energy  T,  integrated  over  energy  transitions  up  to  DT. 
For  DT  > (T  + l)/2,  the  value  of  TCX  is  the  total  cross  section.  In 
addition,  this  subroutine  calculates  certain  partial  cross  sections,  as  follows 

XKAP:  This  is  the  function  which  multiplies  the  stopping  or  cross 

section  number,  to  give  interaction  probabilities,  in  units 
(cm^/g).  That  is,  XKAP  = j (N„/A)  Z/3^,  where  "is 

the  Thomson  cross  section  and  3 = v/c. 

TGE:  This  is  the  glancing  component  due  to  electronic  excitations 

(i.e.,  not  including  vibrations.) 
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TBI: 


This  is  that  part  of  the  glancing  component  for  ionizations 
which  has  as  a factor  the  cut-off  parameter  £n[b(T)]. 


2 

Input  to  this  function  routine  is  the  value  of  BRSQ  = £n(b^R  ),  and. 


through  reference  to  BMOD,  the  correction  term  iin[b^/b(T)] . For  IVB  = 1, 
vibrational  components  are  included  while  for  IVB  = 0 they  are  ignored. 

The  basic  method  involves  use  of  OSCUM  to  determine  the  total  relevant 
oscillator  strength.  To  this  component  is  added  the  knock-on  contribution, 
and,  depending  on  IVB,  the  contribution  due  to  vibrational  excitations. 

Note  that  for  very  small  T,  i.e.,  below  I = 16  eV,  the  value  of  DT 
can  meaningfully  exceed  (T  + I)/2,  which  becomes  inoperable  as  a 

limit.  Note  also  that  the  contributing  glancing  transition  energies  are  not 
allowed  to  exceed  a value  TST,  beyond  which  the  glancing  cross  section  is 
negative. 

Equations  (I.l)  to  (1.8)  in  Appendix  I give  formulae  by  which  the  cross 
sections  are  computed.  We  repeat  below,  using  some  of  the  computer  terms,  the 
important  expressions  and  definitions  for  the  cross  section  functions. 
Quantities  used  in  calculations  of  glancing  collisions  are  as  follows,  with 
U = T/mc^,  I = 1,2,3: 


ARGL  = logK  (2  + U) 


U(U+2) 

(U+1)^ 


(11.18) 


(11.19) 


(11.20) 
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(11.21) 


ASRE(I)  = (1  - nJ 

all  sets  E <DT 
n 

ASRE(3+I)  = X £ (R/E log(l/c  ) (1  - nJ  , (11.22) 

all  sets  E^<DT  n n 

where  the  summations  over  n are  accomplished  by  integration  for  E > 14.2  eV. 
Calculation  of  the  knock-on  contribution  uses  the  following  moments  of  the 
Miller  function: 


(6^  = 16/mc^,  €=  E/mc^,  U = T/mc^,  = 17.5/mc^' 


f. 


DT 


de  6 


I-l 


+ 


1 


2U+1 


U(U+1  )' 


(U+1)' 


. (11.23) 


In  the  case  of  TCX  and  TCXI,  1=1,  while  for  STN,  I = 2,  in 
eqs  (11.19)  to  (11.23). 


6.  TCXI(T,DT,GLI) 

Output  of  this  function  routine  is  the  interaction  probability  at  kinetic 
energy  T for  all  ionizing  processes,  integrated  over  energy  transitions  up 
to  DT.  For  DT  > (T  + I)/2,  the  value  of  TCXI  is  the  total  cross  section 
for  ionization.  In  the  computation  of  this  quantity,  however,  the  glancing 
and  knock-on  processes  are  separately  calculated  and  can  in  principle 'be 
referred  to  separately.  The  argument  GLI  gives  the  glancing  component. 

Input  for  this  function  routine,  and  the  method  of  calculation  are  the 
same  as  for  TCX. 
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7.  STN(T,DT,ISG) 

Output  of  this  function  routine  is  the  stopping  "number"  S/XKAP  for  all 


processes  at  kinetic  energy  T,  integrated  over  energy  transitions  up  to  the 
limit  DT.  For  DT  > (T  + I)/2,  the  value  of  STN  is  the  total  stopping 
number.  The  index  ISG  is  assigned  either  the  value  0 or  the  value  1: 
when  DT  < (T  + I)/2,  both  glancing  and  knock-on  ionizing  transitions  dissipate 
a mean  energy  I = 16  eV  to  atomic  systems  with  each  such  interaction  which 
occurs.  In  a general  way,  this  additional  dissipation  should  be  included  in 
STN  when  DT  is  used  to  identify  "energy  loss";  but  this  dissipation  component 
should  not  be  included  in  STN  whenever  DT  is  used  to  refer  to  (T  + I), 
i.e.,  the  rate  of  energy  loss  due  to  formation  of  secondary  electrons.  In  the 
latter  case  DT  is  then  usually  associated  with  the  maximum  kinetic  energy  of 
secondary  electrons  under  consideration;  hence  electrons  with  greater  kinetic 
energies  are  not  to  be  considered  at  all. 

The  value  of  ISG  is  irrelevant  when  total  stopping  number  is  being 
computed,  i.e.,  for  DT  > (T  + I)/2.  Otherwise  ISG  = 1 signals  inclusion  of 
the  extra  component,  while  ISG  = 0 refers  to  its  exclusion. 

Input  for  this  function  routine  is  the  same  as  in  the  case  of  TCX,  and 
the  method  of  computation  is  essentially  the  same,  except  that  different 
integrals  (identified  by  I = 2)  over  the  transition  energy  are  used  as  indicated 
in  the  notes  for  TCX  given  previously. 

8.  ONCE(TP,DT,XKD,IX) 

Output  of  this  subroutine  is  XKAP"^  times  the  differential  cross 
section  for  various  processes  for  electrons  of  kinetic  energy  TP  undergoing 
energy  loss  DT.  The  output  value  is  XKD.  The  value  of  IX  relates  to  the 
inclusion  or  exclusion  of  secondary  electrons:  IX  = 1 means  that  secondaries 

are  included,  while  for  IX  = 0 or  negative,  secondaries  are  excluded. 
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Input  into  this  subroutine  is  essentially  the  same  as  for  TCX. 

The  calculation  for  glancing  collisions  utilizes  OSCUM,  in  the  mode 
(ID  = 1)  in  which  the  differential  OS  value  for  excitations  or  ionizations 
below  17  eV  is  obtained  by  numerical  differentiation.  The  differential  OS 
values  for  higher  energy  transitions  are  obtained  in  OSCUM  by  interpolation 

on  a predetermined  list.  Cross  sections  for  knock-on  collisions  are  given  by 

2 2 2 
M(j)ller  formula,  with  6^  = 16/mc  , e=  E/mc  , U = T/mc  : 

J_  ^ L 2U+1 

U(U+1)2 

One  should  note  that  the  use  of  a limit  (DT  < TST)  ensures  evaluation 
of  glancing  collision  cross  sections  only  when  they  are  positive.  Also,  for 
DT  > (TP  + I)/2,  glancing  collisions  are  calculated  with  (TP  - DT)  replacing 
DT,  corresponding  to  an  energy  loss  totaling  I plus  the  energy  of  the 
secondary  electron. 


1 + 

€ T+e 


i_ 

m-"j 


11.24' 


(U+1)' 


C.  Subroutines  for  Evaluating  the  Scattering  Integral 
9.  HARDP(N) 

Output  of  this  elementary  tabulation  subroutine  is  a list  HKI(N,J) 
which  is  used  in  the  evaluation  of  the  scattering  integral.  Setting  T = T(N), 
the  n'th  solution  energy,  the  following  values  are  currently  included: 


HKI(N,1 ) 
HKI(N,2) 
HKI(N,3) 


= ARGL  = £n[2(T/R)(2  + T/mc^]  - 3^ 

= 1. 

^ _ (2T/mc^+l ) 

(T/mc^)(l+T/mc^)^ 


HKI(N,4)  = 1/(1  + T/mc^)^  = (1  - 3^) 
HKI(N,5)  = £n(b^/b)  = BMOD 
HKI(N,6)  = £n(l/b^) 

HKI(N,7)  = DTU 

HKI(N,8)  = TB,  as  given  in  FUNCTION  DTU. 


(11.25) 
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These  will  all  be  recognized  as  functions  of  T required  in  cross  section 
formulae.  The  first  six  are  required  sequentially  in  the  order  given  above, 
in  a DO- loop. 

10.  CUMK2(N,NP,A,B,C) 

Output  for  this  subroutine  are  integration  weights  A,B,C  to  be  applied 
to  y(T^p_^  ,Tq)  , respectively.  These  weights 

take  account  of  both  the  energy  region  T^p  to  Tj^  in  the  scattering  integral 
and  the  "scatter-out"  term  K(Tj^)y(Tj^,T^)  usually  placed  on  the  left  side  as  i 
eq  (1).  Note  that  this  region  of  the  scattering  integral  is  complicated  and 
difficult;  the  three  weights  contain  the  contribution  by  the  kernel  function 
k(T',T),  which  has  singular  and  irregular  behavior  for  T'  near  T. 

The  weights  are  obtained  by  Gauss  quadrature,  using  3 moments  and 
3 Dirac  delta  functions  except  when  NP  = 1 , in  which  case  only  2 moments 
and  Dirac  delta  functions  are  used.  The  general  concept  is  described  in 
section  IV,  C in  the  text.  Basical ly, 'one  is  using  eq  (4)  and  making  the 
approximation  = T|.^p  - T^^,  etc.) 

g^E)  »=  o.«(E  - AT|^,p^, ) + B.«(E  - AT,jp)  + yy(E  - AT^p_q  , (11.26) 

such  that 


dE  E'^ 


g,(E) 


a.j  AT 


j 

NP+1 


+ 3.jAT 


j 

NP 


+ 


y.AT 


j 

NP-1 


(11.27) 


One  then  cumulates  all  terms,  with  the  appropriate  factors; 


WT(NP  + 1,N)  = Z + *T^p+l)  ■ 

1 

WT(NP,N)  = Z 3^-f^-(TN  + AT^p)  , (11.28) 

1 

WT(NP-1,N)  = Z Y,-f.(T„  + ATfjp  ,)  . 

i 

We  found  it  desirable  to  include  the  term  K(T|^)y (T|^,T^)  as  a separate 
weight,  WT(N,N),  when  NP  + 1 7^  N.  But  when  NP  + 1 = N,  this  term  was 
incorporated  into  the  integrand  of  eq  (11.27)  by  use  of  the  representation 
K(T^  + E)6(E). 

11.  CUMK3(T,DT,IX) 

Output  for  this  subroutine  is  a set  of  numbers  GL(J),  J = 1,...,6,  on 
which  integration  weights  for  the  scattering  integral  are  based.  These  numbers 
correspond  to  the  following  integrals: 

T+DT  °° 

GL(J)  = J dE  9j(E)  - JdEgj(E)  . (11.29) 

T T 

The  final  term  in  eq  (11.29),  which  is  independent  of  DT  and  is  included 
only  for  glancing  interactions,  always  cancels  because  only  differences 
between  GL  values  appear  in  the  weights.  It  is  used  to  control  round-off 
errors;  note  that  for  large  DT  the  value  of  GL  should  vanish  for  glancing 
interactions.  Further,  DT  can  be  very  large,  with  very  many  integration 
points. 
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The  value  of  IX  determines  whether  or  not  a contribution  due  to  secondary 
electrons  is  to  be  included.  For  IX  = 1 they  are  included,  while  for  IX  = 0 
secondary  electrons  are  excluded.  In  the  case  of  knock-on  interactions,  this 
exclusion  is  effected  by  ensuring  that  for  DT  > (T  + I),  the  upper  limit  in 
the  above  integral  is  .^{2T  +1).  In  the  case  of  glancing  interactions, 
inclusion  or  exclusion  of  secondaries  is  a more  complex  problem  because  the 
cross  section  does  not  explicitly  have  a symmetric  structure.  For  DT  < (T  + I) 
eq  (11.29)  gives  the  relevant  (glancing  interaction)  integrals,  while  for 
DT  > (T  + I),  i.e.,  for  glancing  secondaries,  we  include  a contribution 

.T+DT 

9j(T  + I)  I dE  . (11.30) 

■^2T+I 

One  should  note  however  that  in  all  the  above,  the  limits  of  integration 
can  be  modified  by  the  requirement  that  the  cross  section  never  be  negative,  a 
condition  expressed  through  limits  to  the  value  of  E.  Thus  the  upper  limit  of 
the  integral  in  eq  (11.29)  cannot  exceed  T + DTU,  while  the  lower  limit  of 
the  integral  in  eq  (11.30)  must  be  at  least  TB,  as  evaluated  with  the  DTU 
function  routine. 

12.  WATES(N,IX) 

Output  of  this  subroutine  is  a full  set  of  weights  for  evaluation  of  the 
scattering  integral  at  the  N'th  value  of  T in  the  solution  list.  The 
value  IX  = 1 means  that  these  weights  take  account  of  secondary  electrons, 
while  for  IX  < 0 secondary  electrons  are  excluded. 

Basic  to  this  subroutine  is  the  concept  expressed  in  eq  (9),  the  use  of 
Stieltjes  integration:  values  of  G^.  are  given  in  the  list  6L(I)  which 
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constitutes  the  output  of  CUMK3.  These  values  must  be  turned  into  integration 
weights,  after  which  the  output  of  CUMK2  can  be  included. 

The  procedure  used  (see  section  IV,  C)  is  a bit  indirect.  The  assumption 
is  made  that  the  integrand,  which  is  the  product  f^.  (T  + E)  y(T  + E^T^),  is 
no  more  sharply  varying  than  quadratic  in  an  arbitrary  variable  u such  that 
u = - 1,0,1  at  three  successive  points  T^,  , and  Tj+2'  these  same 

three  points  the  basic  integration  variables  G^.  have  values  G^.  (1), 
and  G^.  (3).  We  then  write 


(IT. 31) 


We  next  fit  G^.  (u)  by  a quadratic  function  of  u,  from  which  (dG^/du)  can 
be  approximately  evaluated; 

2 

G^-  (u ) = a + bu  + cu  , 

G.  (-1 ) = a - b + c , 

etc.  (11.32) 


Three-point  integration  weights  can  then  be  applied  for  evaluation  of  the 
integral.  While  this  integration  procedure  may  not  appear  so,  it  is  exact 
if  the  basic  assumption  is  obeyed,  that  f^.y  can  be  exactly  represented  by 
a quadratic  function  of  the  variable  u. 

As  in  the  case  of  CUMK2,  the  cumulation  over  i values  with  the  f^. 
factors,  and  combination  of  different  weights  for  the  same  point  must  be 
carried  out. 
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We  should  note  that  general  3-point  or  Simpson  rule  calculations  are 


suitable  to  an  odd  number  of  points,  whereas  the  scattering  integral  will 
involve  odd  and  even  numbers  of  points  alternately.  For  cases  in  which  an 
even  number  of  points  occurs,  the  interval  from  TY(1)  to  TY(2)  is  separ- 
ately evaluated  as  follows:  The  procedure  previously  described  is  followed 

for  points  TY(1),  TY(2),  and  TY(3),  except  that  the  integration  over  the 
variable  u is  taken  only  from  -1  to  0.  The  resulting  three  weights  are 
included  with  the  same  cumulations  and  factors  already  described. 

The  subroutine  closes  by  calling  CUMK2  and  including  weights  for  the 
final  integral,  including  the  "scatter-out"  term. 

13.  C0RRX(X,UI2,UI1 ,UI0,A,B,C) 

Output  of  this  small  correction  subroutine  are  three  values  A,  B,  and 
C,  which  are  essentially  differences  between  B-poi'nt  integration  weights  based 
on  the  assumption  of  a continuous  integrand,  and  Gauss  quadrature  weights  which 
can  take  account  of  a discontinuity  in  the  same  interval. 

The  three  points  of  the  integration  interval  are  UI2,  UIl,  and  UIO, 

in  that  order.  The  value  of  X falls  somewhere  in  the  interval,  this  being 

the  location  of  the  discontinuity,  with  X = 0 corresponding  to  UI2.  Size 

_2 

of  the  discontinuity  is  I , as  shown  in  the  sketch  below: 


The  function  sketched  is  a constant  times  a unit  function. 
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Such  a discontinuity  occurs  in  the  differential  knock-on  cross  section; 

and  we  take  account  only  of  the  discontinuous  part  of  the  largest  term, 

7 7 

E = I , evaluated  at  the  discontinuity.  Integration  as  in  CUMK3  here 
involves  the  variable  (see  eq.  (4)) 

dG  = 1 dX  , (11.33) 

and 


The  procedure  used  is  that  of  constructing  weights  at  the  three  points  by 
straightforward  application  of  three  point  numerical  integration  as  discussed 
in  WAXES.  But  we  also  calculate  weights  by  applying  Gauss  quadrature  (in  the 
variable  X)  to  the  step  function  by  Dirac  delta  functions  located  at  the 
three  points  indicated. 

Differences  between  the  two  sets  of  weights  give  the  values  of  A,  B, 
and  C. 

The  reason  for  this  calculation  is  our  impression  that  the  Gauss  quadra- 
ture approach  is  significantly  more  accurate  for  any  interval  containing  the 
discontinuity  generated  at  the  17.5  eV  cut-off  of  the  knock-on  transitions. 

D.  General  Evaluation  of  the  Slowinq-Down  Spectra 

14.  TLIST(NPH) 

Output  of  this  subroutine  is  a solution  list  of  energies,  N3  in  number, 
which  are  stored  in  the  array  TY(I).  A peculiarity  of  this  list  is  that  for 
T > T^/2,  the  quantity  TY ( I ) corresponds  to  (T^  - T),  while  for  T < Tq/2, 
TY(I)  = T.  The  reason  for  this  distinction  is  mostly  the  avoidance  of  round- 
off error  in  the  computations. 
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The  Input  constant  NPH  gives  the  factor  by  which  T changes  (in  the 


intermediate  energy  region)  between  successive  values,  viz 

T - t = <To  - Vl)  • • (”-35) 

The  input  constant  N1  determines  how  many  finely  spaced  intervals  are 
utilized.  These  finely  spaced  intervals  are  separated  approximately  as  the 
basic  tabulations  of  OSCUM,  and  as  given  in  an  input  list  EOS.  But  the 
TY(I)  are  located  between  the  EOS  values  according  to  the  formula 

TY(I)  = {EOS(I)  • E0S(I  + 1)}^^  . (11.36) 

After  N1  of  these  values,  two  buffer  intervals  of  10  eV  are  inserted 
separating  the  fine  spacing  from  the  coarse  spacing.  The  Gauss  quadrature  is 
arranged  to  extend  over  an  energy  region  exceed i ng  11.2  eV.  When  (T^  -j  - T^) 
corresponds  to  the  first  of  these  buffer  intervals,  fine  spacing  can  resolve 
the  peak  region  of  the  oscillator  strength  distribution;  while  for  (T^  .j  - T^^) 
corresponding  to  the  second  of  these  intervals,  the  peak  region  is  covered  by 
the  Gauss  quadrature  region,  now  20  eV  wide. 

Coarse  spacing  begins  with  an  interval  of  20  eV  and  increases  at  the 
rate  already  given. 


AT  ^t/AT  = 2 
n+1  n 


1/NPH 


(11.37) 


Intervals  increase  in  spacing  until  exceeded,  at  some  value  T^. 

For  this  value,  one  writes  N2  = n,  and  TY(N2)  = Tq/2.  The  value 
TY(N2  - 1)  is  then  readjusted  to  give  a more  even  sequence  near  Tq/2. 
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For  n > N2,  the  values  of  TY(n)  repeat  in  reverse  order  those  which 


precede,  according  to  the  rule 


TY(n)  = TY[2-N2-n]  . (11.38) 

This  has  the  consequence  that  the  intervals  decrease  in  width  progressively 
towards  low  electron  kinetic  energy. 

Finally,  intervals  may  be  added  as  necessary.  In  particular,  10  to  20 
additional  values  are  usually  added  below  the  lowest  electronic  transition 
energy  to  cover  the  region  in  which  vibrational  transitions  are  dominant. 

(No  attempt  is  made  in  the  interval  near  to  perform  corresponding 
calculations. ) 

One  should  note  that  for  low  source  energies  T^,  the  sequence  just 
described  is  interrupted  whenever  the  tabulated  value  recorded  in  TY  exceeds 
T^/2.  Below  T^/2,  however,  the  procedure  is  as  already  outlined. 

Both  N2  and  N3  are  output  constants  that  are  recorded  in  memory 
locations  available  to  other  subroutines. 

15.  CXPRT(KX) 

The  name  for  this  subroutine  is  short  for  "cross  section  printout." 
Output  determined  by  this  routine  includes  the  stopping  number,  the  total 
cross  section  and  various  components,  including  the  total  ionization  cross 
section,  and  the  function  XKAP  which  multiplies  stopping  number  and^total 
cross  section  "number"  to  give  stopping  powers  and  interaction  probabilities 
per  unit  path  length. 

For  KX  = 1 , the  differential  cross  section  is  printed  as  wen, 
normalized  to  be  compared  with  a Platzman  plot.  For  KX  = 0 (or  negative), 
only  the  integral  cross  sections  are  recorded. 
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16.  SOLVE(IS) 

Output  of  this  subroutine  is  a solution  list  for  the  index  IS;  the 


values  are  recorded  in  the  array  S0L(I,2-IS).  Interpretation  of  the 
different  possible  cases  is  as  follows: 

IS  = 1 : The  total  slowing-down  spectrum  including  secondaries  of  all 

orders; 

IS  = 0:  primary  electrons  only 

IS  = -1:  1st  generation  secondaries  only; 

IS  = -J:  J'th  generation  secondaries  only. 

For  IS  = 1,  the  integration  weights  required  are  those  based  on 
inclusion  of  secondary  electrons.  For  IS  = 0,  the  integration  weights 
are  zero  for  secondary  electron  components.  For  negative  IS  values,  the 
integration  weights  which  are  appropriate  are  di fferences  between  the  above, 
i.e.,  weights  only  for  secondary  electrons,  which  are  used  in  tandem  with 
weights  only  for  primary  electrons. 

For  negative  IS  values,  the  two  required  sets  of  weights  are  put  above 
and  below  the  diagonal  element  of  the  WT  array,  with  the  diagonal  reserved 
for  WT(N,N). 

The  calculation  proceeds  by  setting  up  an  array  of  integration  weights  as 
called  for.  However,  the  different  generations  of  secondary  electrons  are 
calculated  in  order;  and  all  use  the  same  array  of  weights,  which  is  thus 
determined  only  once,  and  used  for  all  J < - 1. 
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Next,  the  calculation  proceeds  by  evaluating  the  source  contribution 
(ONCE),  the  scattering  integral  contribution,  and  finally,  by  use  of  WT(N,N), 
the  value  of  y.  Because  the  different  generations  of  secondaries  are  computed 
in  chain  fashion,  the  solution  for  the  preceding  generation  is  retained  in  a 
special  list,  y(I),  for  further  reference.  (This  is  convenient  rather  than 
essential,  because  this  information  is  also  in  the  general  solution  array  SOL.) 


17.  ECHEK(N,IX) 

Output  is  the  fraction  of  the  input  energy  which  is  dissipated  in 
the  form  of  electrons  of  kinetic  energy  lower  than  the  N'th  solution  energy 
T,  plus  the  energy  deposited  with  molecules,  due  both  to  excitations  and 
ionizations.  Also  computed  is  the  number  of  electrons  making  the  transition 
from  kinetic  energy  above  T|^  to  kinetic  energy  below 

Because  all  energy  is  dissipated  regardless  of  the  criterion  used,  energy 
conservation  requires  a value  of  unity  for  the  fraction  mentioned  above. 
Similarly,  number  conservation  requires  a value  of  unity  for  the  number  of 
electrons  crossing  the  cut-off  for  T such  that  > (T  - I)/2  or  for 

calculations  based  on  the  spectrum  of  primaries  only. 

The  index  IX  merely  selects  the  spectrum  whose  energy  fraction  is  to  be 
computed.  Note  that  energy  and  number  dissipation  can  be  computed  not  only 
for  the  total  slowing-down  spectrum,  but  also  for  its  various  components, 
viz.,  primaries  only,  or  n'th  generation  secondaries. 

Let  us  write  k(T',6)  for  ^ ^ ' where  e=  T'  - T,  where  k is  a 

function  earlier  called  XKAP, 


k{1' 


3^(T 


(11.39) 
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We  then  define  various  derived  quantities  as  follows: 


JC(T',x)  = f d6  k(T',6),  i.e.,  TCX(T' ,x,IVB) 

-'o 


K(r  ) =Jt(T' ,(T'+I)/2),  also  TCX(TMMVB) 


^(T',x)  = J de  ek(T',e),  i.e.,  STN(T',x,0) 


S(T')  =^(T',(T'+I)/2),  also  STN(T',T',0)  , 


J*(T',x)  =^(r,x)  + I[K(T')  -J((T',x)],  i.e.,  STN(T',x,l) 


(11.40) 


(11.41) 


(11.42) 


(11.43) 


(11.44) 


Then  eq  (7)  is  interpreted  numerically  with  the  following  identifications: 
ydM^)  =ay(T',T^)/K(T')  , (11.45) 

S(T',T)  = k(T')  _/*(!', T+I),  T’>2T+I  , 

(11.46) 

= k(T')  { + T'[K(T')  -JC(T' ,T'-T)]},  T-  < 2T  + I . 
Similarly,  eq  (8)  can  be  interpreted  as  follows: 

Kp(T',T)  = k(T')[K(T')  - J{(T' ,T'-T)],  T'  < 2T  + I and  T'  < , (11.47) 

K3(T',T)  = k(T')[K(T')  - J((T',T+D],  T'  > 2T  + I . (11.48) 
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In  both  cases,  solution  is  for  0/(T',T^).  But  to  check  number  conservation, 
it  is  more  useful  to  take  = 0. 

The  integrands  of  eqs.  (7)  and  (8)  become  very  large  for  T'  T.  In 
principle  this  suggests  an  application  of  Gauss  quadrature.  But  we  prefer 
here  simply  to  add  and  subtract  integral  quantities  with  two  properties: 
a)  the  modified  integrands  tend  to  vanish  as  T'  ->  T;  and  b)  the  subtracted 
integral  term  can  be  accurately  and  readily  evaluated.  Our  choice  of  supple- 
mentary term  for  evaluation  of  the  energy  conservation  equation  is  the 
following  identity: 


Correspondingly , we  use  the  following  term  in  the  number  conservation 
equation. 


(3T+I)/2 


dT-  yd,!^)  T[K(T)  -.X(T,T'-T)]  = T y(T,T^)  S(T)  . (11.49) 


T 


(3T+I)/2 


dT'  y(T,T^)  [K(T)  -jf(T,T'-T)]  = y(T,T^)  S(T)  . (11.50) 


T 


Equations  (7)  and  (8)  then  take  the  forms 


T 

I 


''  0 


T y{T,iy  S(T)  = Tg 


dT'  y(T',Ty  k(T')  S(T',T) 


- , (11.51) 


T 


y(T,Tjj)  S{T)  = 1 +1  dT'  y(T',Tg)  k(T')  Kj(T',T) 


2T+I 

2T+I 


(11.52) 
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If  one  combines  terms  in  braces  into  single  integrals,  one  finds  that  the 
integrands  then  approach  zero  in  the  limit  of  T'  T.  One  remaining  question 
pertains  to  whether  the  approach  to  zero  may  not  be  so  slow  that  the  contri- 
bution of  the  final  interval  is,  say,  underestimated;  this  can  in  principle  be  a 
significant  problem  because  of  an  often  sizable  contribution  due  to  this 
interval.  To  keep  track  of  this,  an  extrapolation  of  the  integrand  to  a 
finite  value,  rather  than  zero,  has  been  made,  with  an  associated  estimate  of 
the  corresponding  modification  to  the  integrated  result.  On  the  basis  of  this 
estimate,  the  problem  did  not  appear  to  be  serious. 

The  above  addition-subtraction  is  not  carried  out  below  the  lowest 
electronic  transition  energy,  taken  to  be  THR  = 11.25  eV.  In  that  low 
energy  region,  evaluation  of  the  integrands  is  carried  out  directly. 

The  outline  of  the  computation  is  as  follows:  First,  a list  of  integration 

weights  is  computed,  using  the  Stieltjes  integration  method  described  in 
section  II. 11.  This  method  was  not  essential;  but  the  intervals  vary 

in  size,  and  this  is  a suitable  way  to  take  account  of  these  variations. 

Next,  the  numerical  integration  sum  is  evaluated;  and  evaluations  of  the 
addition-subtraction  and  the  source  term  follow.  These  three  terms,  their 
sum,  and  the  extrapolation  fractional  increment  used  to  check  the  addition- 
subtraction  procedure  are  recorded  for  possible  output.  A sixth  quantity  is 
also  evaluated,  namely  the  number  of  electrons  crossing  the  specified  energy 
cut-off  T|,^. 

In  addition,  three  other  quantities  are  routinely  evaluated.  One  of 
these  is  the  value  of  1/W(T),  number  of  ions  per  eV  generated  by 
y(T'),  r > T.  The  second  is  the  total  yield  of  excitations.  The  third  is 
the  factor  by  which  the  solution  niust  be  modified  to  bring  about 

accurate  energy  conservation.  These  last  three  quantities  are  stored  separately 
from  the  other  6 quantities. 
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18.  ELSPEL 


Output  of  ELSPEL  constitute  the  main  results  of  the  computation, 
including  stopping  powers  and  cross  sections,  the  value  for  W,  the  yield  of 
excitations  of  all  types,  the  slowing-down  spectrum  and  its  components,  and 
the  results  of  the  energy-and  number-conservation  tests. 

Input  quantities  are  as  follows: 


NEX:  the  number  of  source  energies, 

IX:  the  number  of  solution  lists.  IX  = 1 here  corresponds  to 

the  slowing-down  spectrum  alone.  IX  = 2 calls  for  a 
calculation  also  of  the  primaries  as  a separate  component. 

IX  = 3 calls  for  the  1 st-generation  secondaries  as  well,  etc. 


SORS(I):  the  list  of  source  energies,  in  eV. 

BCOR(I):  see  BMOD 

NO,  N1 : see  TLIST. 

B:  the  value  of  b . 


In  addition,  NPH  determines  the  number  of  solution  energies  and  their 
separation,  as  discussed  in  TLIST. 

The  main  DO-loop  permits  the  computations  to  be  repeated  for  different 
NPH  values  (see  also  TLIST),  to  check  effects  due  to  size  of  integration 
mesh.  In  this  DO-loop  the  solution  energies  are  first  fixed,  by  calling 
TLIST.  Then  CXPRT  is  called  to  record  various  stopping  powers  and  cross 
sections.  Next  HARDP  performs  the  preliminary  computations  of  necessary 
energy-dependent  parameters.  SOLVE  is  then  called  to  compute  the  slowing- 
down  spectrum;  and  this  computation  is  followed  by  print-out  of  results. 
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together  with  related  quantities  including  the  integrand  for  computation  of 
(1/W),  Finally,  ECHEK  performs  the  energy-  and  number-conservation  tests, 
and  values  for  W and  any  for  excitation  yields. 


The  whole  computation  can  then  be  repeated  at  other  value  of 
NEX  >1,  or  possibly  with  a different  set  of  integration  points. 
Other  special  cross  section  routines  are  sometimes  inserted  in 


ELSPEL. 
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Table  1 
Table  2 

Table  3 
Table  4 

Table  5 
Table  6 


Table  Captions 


Integration  energies  for  a calculation  with  T^  = 1 keV. 

Values  for  energy  and  number  conservation  integrals  above 
different  T,  expressed  as  fractions  of  the  source  energy  T^, 
for  T^  = 30  keV.  For  the  column  on  the  far  right,  integra- 
tion was  carried  out  only  over  the  primary  component. 

Ratios  y (T,T^)/y (T,l  MeV),  for  various  T,T^. 

Experimental  and  theoretical  values  for  W,  FE,  F,  and  DS 
refer  to  Fowler  equation,  (Fano)  energy  balance,  and  degradation 
spectrum  methods  of  calculation  [55]. 

Volts  per  ion  pair  (W)  as  computed  for  different  source  electron 
kinetic  energies. 

Some  comparisons  with  data  from  reference  [15]. 
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Table  1.  Integration  energies  for  a calculation  with  T = 1 keV. 
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29.42  53  404.4  80  31.74  107  11.43 

30.56  54  500  81  30.56  108  10.93 
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Table  3 

. Ratios  y(T,T^ 

)/y(T,l  MeV), 

for  various 

T’To- 

To,keV 

T = 9.9 

12.3 

17.3 

22.6 

.10 

1.376 

1.470 

1 .813 

2.197 

.30 

1.132 

1.163 

1.279 

1 .401 

.50 

1.112 

1.134 

1.211 

1 .288 

1.0 

1.0775 

1 .0922 

1 .144 

1.190 

3.0 

1.0564 

1.0653 

1 .0965 

1 .122 

5.0 

1.0490 

1.0564 

1 .0824 

1 .103 

10. 

1.0394 

1.0453 

1.0658 

1.0822 

30. 

1.0269 

1.0308 

1 . 0448 

1 .0558 

100 

1 .0159 

1 .0181 

1.0265 

1 .0304 

300 

1 . 0080 

1.0090 

1 .0133 

1 .0166 

67 


Table  4.  Experimental  and  theoretical  values  for  W (see  also  Whyte  [56]  and  Booz  and  Ebert  [21]). 
FE,  F,  and  DS  refer  to  Fowler  equation,  (Fano)  energy  balance, 
and  degradation  spectrum  methods  of  calculation  [55]. 
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Table  5.  Volts  per  ion  pair  (W)  as  computed  for 
different  source  electron  kinetic  energies 


T,  keV 

w 

.024 

69.4 

.027 

73.2 

.030 

67.4 

.035 

56.5 

.05 

50.7 

.07 

47.1 

.1 

44.2 

.3 

39.27 

.5 

37.44 

1 . 

36.24 

2. 

36.03 

3. 

36.05 

5. 

36.08 

10. 

36.20 

30. 

36.40 

100. 

36.59 

300. 

36.75 

1000. 

36.90 
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Table  6.  Some  comparisons  with  data  from  reference  [1]. 


T(eV) 


y (cm/eV) 


Okazaki -Sato 


10 


,26  X 10 


-3 


.46  X 10' 


10^ 


23  X 10‘ 


.43  X 10' 


10^ 


80  X 10 


-3 


1.22  X 10‘ 


Table  Captions  - Appendix  I 


Table  I. 
Table  I. 
Table  I. 
Table  I. 

Table  I. 
Table  I. 


Table  I 


1 Oscillator  strengths  for  discrete  excitations  [3],  based  on  Allison 
and  Dalgarno  [59]  and  Namioka  [29]. 

2 Tabulations  of  moments  of  the  OSD,  as  defined  in  eqs  (II. 1)  and 
(II. 2). 

, ""i 

3 A tabulation  of  low  energy  correction  functions  for  ionization 
(£n  b^/b)  using  b^  = 3,  and  excitation  (£n  1/b.j). 

4 Comparisons  of  total  ionization  and  excitation  cross  sections  for 
electrons  with  kinetic  energies  below  1000  eV.  See  references 
[3,33].  (The  number  of  digits  given  is  for  comparisons;  neither 
type  of  cross  section  is  known  to  any  such  accuracy.) 

5 Vibrational  cross  sections  for  the  two  lowest  excitations,  used  in 
conjunction  with  eqs  (8)  and  (9)  and  due  to  Trajmar,  et  al  [64]. 

6 Values  of  <I),  the  mean  excitation  potential  for  molecular  hydrogen. 


Table  Captions  - Appendix  II 


I.l  Description  of  ELSPEL  subprograms.  Abbreviations  are  as  follows: 
c.s. — cross  section,  c . c . s . --combined  cross  section,  sc.int-- 
scattering  integral,  eval . --evaluation 


% 
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Table  I.l  Oscillator  Strengths  for  discrete  excitations  [3], 
based  on  Allison  and  Dalgarno  [59],  and  Namioka  [29]. 


E„.  eV 

Lyman 

f 

n 

E , eV 
n 

Werner 

f 

n 

11.213 

.001709 

12.317 

.05246 

11.377 

.005859 

12.605 

.08027 

11 .535 

.01170 

12.875 

.07694 

11.674 

.01776 

13.129 

.06030 

11.846 

.02277 

13.367 

.04269 

1 1 . 982 

.02602 

13.589 

.02863 

12.132 

.02736 

13.794 

.01873 

12.268 

.02705 

13.984 

.01211 

12.405 

.02553 

14.156 

.007822 

12.524 

.02326 

14.319 

.005060 

12.679 

.02059 

14.446 

.003280 

12.790 

.01785 

14.560 

.002104 

12.902 

.01522 

14.650 

.001290 

13.026 

.01281 

14.710 

.000616 

13.139 

.01068 

13.238 

.008835 

13.305 

.007271 

13.450 

.005962 

13.549 

. 004878 

13.648 

.003986 

13.735 

.003258 

13.822 

.002664 

13.909 

.002180 

13.996 

.001787 

14.071 

.001467 

14.150 

.001207 

14.220 

.0009933 

Total 

.3107 

Total 

.J923 
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Table 


.3.  A tabulation  of 

low  energy  correction  functions  for  ionization 

(An  b^/b)  using 

b =3,  and  excitation 

(£n  1/b^ ) . 

T,  eV 

£n  (b  /b) 

' oo'  ' 

£n  (1/b-j  ) 

1414 

.00 

.000 

1000 

.13 

.000 

707 

. .25 

.000 

500 

.40 

.000 

353.6 

.62 

.000 

250 

.90 

.010 

176.8 

1.23 

.175 

125 

1.54 

.408 

88.4 

1.78 

.614 

62.5 

1.93 

.819 

44.2 

1.99 

1.037 

31.25 

2.02 

1.185 

22.1 

1.92 

1 .210 

18.58 

1.79 

15.625 

2.40 

1 .100 

14.00 

.98 

11.05 

3.20 

.81 
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Table  1.4.  Comparisons  of  total  ionization  and  excitation  cross  sections 
for  electrons  with  kinetic  energies  below  1000  eV.  See  refs.  [3,33]. 
(The  number  of  digits  given  is  for  comparisons; 
neither  type  of  cross  section  is  known  to  any  such  accuracy.) 
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073 


Table  1.5,  Vibrational  cross  sections  for  the  two  lowest  excitations, 
used  in  conjunction  wi  h eqs.  (8)  and  (9)  and  due  to  Trajmar,  et  al  [64], 


Ta 

4rra^R 

0 


T,  eV 

0.54  eV 

1.08  eV 

1.3 

0 

0 

1.6 

.0081 

.000110 

2.0 

.0162 

.000560 

2.5 

.0251 

.00147 

3.0 

.0321 

.00220 

3.5 

.0374 

.00276 

4.0 

.0389 

.00324 

5.0 

.0383 

.00360 

6.0 

.0357 

.00382 

7.0 

.0329 

.00391 

8.0 

.0254 

.00391 

9.0 

.0229 

.00320 

10.0 

.0181 

.00205 

14.0 

.00819 

.000737 

20.0 

.00585 

.000199 
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Table  1.6  Values  of  <I),  the  mean  excitation  potential  for  molecular  hydrogen. 
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Table  II. 1.  Description  of  ELSPEL  subprograms.  Abbreviations  are  as 
follows:  c.s. --cross  section,  c . c . s . --combi ned  cross 

section,  sc. int. --scattering  integral,  eval . --eval uation 


SUBPROGRAM  TYPE  DEPENDENCE  PURPOSE  OR  CONTEXT 


OSCUM 

c.s. 

VIB 

c.s. 

BMOD 

c.s. 

DTU 

c.s. 

BMOD 

TCX 

c.c.s. 

all  c.s. 

TCXI 

c.c.s. 

all  c.s. 

STN 

c.c.s. 

all  c.s. 

ONCE 

c.c.s. 

all  c.s. 

HARDP 

sc. int. 

DTU 

CUMK2 

sc . int. 

all  c. s. , 
HARDP, TCX 

CUMK3 

sc. int. 

all  c.s. 

CORRX 

sc. int. 

WATES 

sc .int. 

all  other 
sc. int. 

TLIST 

eval . 

ECHEK 

eval . 

all  c.c.s. 

SOLVE 

eval . 

WATES, ONCE, 
ECHEK 

CXPRT 

eval . 

all  c.c.s. 

ELSPEL 

master 

al 1 eval . 

cumulative  oscillator  strength 
vibrational  cross  sections 
low  energy  glancing  c.s.  parameters 
glancing  c.s.  limits 

cumulative  cross  sections 
cumulative  ioniz.  c.s. 
cumulative  stopping  power 
differential  cross  sections 

factors,  parameters  for  integration 
Gauss  quadrature,  final  interval 

cumulative  integrations 
correction  for  cut-off  discontinuity 
array  of  integration  weights 


prepare  list  of  solution  energies 
evaluate  energy,  number  conservation 
evaluate  slowing-down  spectrum 

evaluate  total  c.s.,  stopping  powers 

direct  all  computations,  output 
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Figure  Captions 


Fig.  1 
Fig.  2 

Fig.  3 

Fig.  4 
Fig.  5 

Fig.  6 
Fig.  7 

Fig.  8 
Fig.  9 

Fig.  1 
Fig.  1 


Differential  cross  sections  for  electrons  of  kinetic  energy  T as 
a function  of  energy  loss,  T -T.  Secondary  electrons  comprise  that 
part  of  the  spectra  below  (T°-16)/2  eV. 

Differential  cross  sections  for  secondary  electrons  of  final  energy 
E due  to  500  eV  and  100  eV  electrons  in  H . Solid  curves  are 
calculated  from  our  cross  section  representations  (R  = 13.6  eV). 

The  500  eV  curve  agrees  closely  with  Gerhart's  data  from  his  Fig.  7, 
for  E ^ > .05.  Circles  are  Opal,  Beaty,  and  Peterson  data  [30] 
divided  by  1.65.  The  + values  are  unpublished  data  due  to  Dubois 
and  Rudd  [31]. 

The  stopping  power  for  H^.  .The  dashed  line  gives  results  for  the 
Bethe  stopping  power  formula  using  19.22  eV  for  the  mean  excitation 
potential.  The  circles  give  recent  data  due  to  Green  [32]. 

Total  cross  section  for  inelastic  interactions  of  fast  electrons  in 
H^.  Also  given  is  the  fraction  due  to  ionizing  interactions. 

Relative  contributions  to  the  energy  dissipated  to  kinetic  energies 
below  16  eV,  as  a function  of  the  kinetic  energy  of  the  source 
electrons. 

The  Lewis  effect  for  a 1 MeV  source.  S(T)  is  the  stopping  power. 
Solution  intervals  in  the  Lewis  effect  region  varied  from  ~ 0.4  eV 
at  (Tq-T)  ~ 11  eV,  to  ~ 2 eV  at  (T  -T)  ~ 60  eV.  The  solution  was 
modified  to  improve  energy  conservation  in  the  region  between  the 
arrows  at  the  bottom. 

Suitably  scaled  electron  slowing-down  spectra  for  source  energies  of 
30  eV,  300  eV,  3 keV,  and  30  keV.  S(T)  is  the  stopping  power.  For 
the  30  keV  source,  primaries  and  1st  generation  secondaries  are 
indicated  with  dashed  lines.  The  ratio  k./S  is  also  given  (heavy 
dashed  line),  as  are  several  relevant  reciprocal  stopping  numbers 
(at  the  right).  Below  about  12  eV,  data  for  the  three  higher 
source  energies  are  so  close  together  that  they  are  all  represented 
by  a single  curve. 

Volts  per  ion  pair  (W)  as  a function  of  electron  kinetic  energy  in 
H,, . Comparison  points  are  given  for  Combecher's  experimental  values 
[57],  and  calculated  values  from  references  [2]  and  [3]. 

Volts  per  ion  pair  (W)  in  for  electron  kinetic  energies  greater 
than  1 keV.  The  scale  is  enlarged  so  that  small  differences  of 
trend  can  be  readily  seen.  At  the  right,  a change  of  ]%  is  exhibited 
The  circles  give  values  due  to  Garvey,  Porter  and  Green  [2],  as  read 
from  their  curve.  Green,  Jackman,  and  Garvey  (see  [2])  give  a few 
numerical  values  which  are  in  agreement  to  three  significant  figures. 

Comparison  of  values  for  G = 100/W.  The  solid  line  represents  our 
calculations.  Circles  give  Combecher's  experimental  data  [57]. 

Recent  calculations  due  to  Douthat  are  also  shown  [54]. 

1 The  slowing-down  spectrum  of  10  keV  electrons  in  H^,  compared  with 
values  due  to  Douthat  [1,54]  and  due  to  Garvey,  Porter,  and  Green  [2] 
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Fig.  2 Differential  cross  sections  for  secondary  electrons  of  final  energy 
E due  to  500  eV  and  100  eV  electrons  in  . Solid  curves  are 
calculated  from  our  cross  section  representations  (R  = 13.6  eV). 

The  500  eV  curve  agrees  closely  with  Gerhart's  data  from  his  Fig.  7, 
for  E~^  ^ .05.  Circles  are  Opal,  Beaty,  and  Peterson  data  [30] 
divided  by  1.65.  The  + values  are  unpublished  data  due  to  Dubois 
and  Rudd  [31]. 
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The  stopping  power  for  H . The  dashed  line  gives  results  for  the 
Bethe  stopping  power  formula  using  19.22  eV  for  the  mean  excitation 
potential.  The  circles  give  recent  data  due  to  Green  [32]. 
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Total  cross  section  for  inelastic  interactions  of  fast  electrons  i 
Also  given  is  the  fraction  due  to  ionizing  interactions. 
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Fig.  5 Relative  contributions  to  the  energy  dissipated  to  kinetic  energies  below  16  eV,  as 
function  of  the  kinetic  energy  of  the  source  electrons. 
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Fig.  8 Volts  per  ion  pair  (U)  as  a function  of  electron  kinetic  energy  in 
H . Comparison  points  are  given  for  Combecher's  experimental  values 
[57],  and  calculated  values  from  references  [2]  and  [3]. 
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Fig.  10  Col.iparison  of  values  for  G = 100/W.  The  solid  line  represents  our 
calculations.  Circles  give  Combecher's  experimental  data  [57]. 
Recent  calculations  due  to  Douthat  are  also  shown  [54]. 


, cm/eV 


Fig.  11  The  si owi ng-down  spectrum  of  10  keV  electrons  in  H , compared  with 

values  due  to  Douthat  [1,54]  and  due  to  Garvey,  Porter,  and  Green  [2]. 
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Figure  Captions  - Appendix  I 


Fig.  I.l  The  ionization  efficiency  function  as  reproduced  from  [3].  The 

solid  curve  is  due  to  Berkowitz  [61].  The  dashed  line  is  a smoothed 
form  due  to  Gerhart.  The  dash-dot  line  is  due  to  eg  (1.3). 

Fig.  1.2  Gerhart's  $(E)  (heavy  curve)  is  reproduced  from  his  Fig.  6 [3]. 

The  light  curve  expresses  a unit  function  with  cut-off  at  17.5  eV. 
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The  ionization  efficiency  function  as  reproduced  from  [3].  The  solid  curve  is  due  to 
Berkowitz  [61].  The  dashed  line  is  a smoothed  form  due  to  Gerhart.  The  dash-dot  line  is 
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Gerhart's  4>(E)  (heavy  curve)  is  reproduced  from  his  Fig.  6 [3].  The  light  curve  expresses 
a unit  function  with  cut-off  at  17.5  eV. 
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